Author

Qian Hui TAN

Published

October 5, 2023

1 02_DE

In this notebook, we run quick DE volcano plots, p-value histograms and pairwise cross plots for our conditions.

Code
source("../../scripts/R/de_functions.R")
source("../../scripts/R/go_kegg_functions.R")

suppressPackageStartupMessages({
  library(DESeq2)
  library(GenomicFeatures)
  library(RColorBrewer)
  library(gplots)
  library(biomaRt)
  library(grid)
  library(gridExtra)
  library(ggplot2)
  library(lattice)
  library(reshape)
  library(geneplotter)
  library(ggrepel)
  library(limma)
  library(tidyverse)
  library(eulerr)
  library(gghighlight)
  library(clusterProfiler) # install from github
  library(org.Dm.eg.db)
})

1.1 Functions

Code
## If the output directory doesn't exist, create it
if(!dir.exists("../output/02_DE")){
  dir.create("../output/02_DE")
}

output_dir <- "../output/02_DE/"

Loading data:

Code
dds = readRDS("../output/01_QC/dds_filt.RDS")
ensembl.genes= readRDS("../output/01_QC/ensembl_genes.RDS")
lfc.threshold = log2(1.5)


dim(dds)
[1] 15655    20

Remove lowly expressed genes:

Code
# Set control as reference
#colData(dds)$condition <- relevel(colData(dds)$condition, ref = "WT")


# Removing lowly expressed genes, only to be done once at the start of the differential expression step
filter = apply(counts(dds, normalized = TRUE), 1, function(x){ mean(x) >= 10 })
dds = dds[filter, ]
dim(dds)
[1] 11145    20

Run DESeq:

Code
dds <- DESeq(dds, test = "Wald", parallel = TRUE)
using pre-existing size factors
estimating dispersions
gene-wise dispersion estimates: 6 workers
mean-dispersion relationship
final dispersion estimates, fitting model and testing: 6 workers
Code
rld <- rlog(dds, blind = FALSE) # blind = FALSE for downstream DE

design(dds) 
~condition

List of samples in here

Code
colData(dds)$sample_id
 [1] "Ctrl_1"    "Ctrl_2"    "W4_OE_1"   "W4_OE_2"   "O1_OE_1"   "O1_OE_2"  
 [7] "O1W4_OE_1" "O1W4_OE_2" "O2_KO_1"   "O2_KO_2"   "W4_KO_1"   "W4_KO_2"  
[13] "O1_KO_1"   "O1_KO_2"   "O1W4_KO_1" "O1W4_KO_2" "O2W4_KO_1" "O2W4_KO_2"
[19] "O1O2_KO_1" "O1O2_KO_2"

1.2 RLE plot

Code
assay(rld) %>%
  as.data.frame() %>%
  gather() %>%
  ggplot(aes(x=key, y=value)) + 
  geom_boxplot() + ylab("Log Expression ") + xlab("Sample") +
  theme(axis.text.x = element_text(angle = 90))

Run Wald test

Code
wald_dds = nbinomWaldTest(dds)
found results columns, replacing these

Check names of coefficients

Code
res_names <- resultsNames(wald_dds)

res_names
 [1] "Intercept"                 "condition_O1_KO_vs_Ctrl"  
 [3] "condition_O1_OE_vs_Ctrl"   "condition_O1O2_KO_vs_Ctrl"
 [5] "condition_O1W4_KO_vs_Ctrl" "condition_O1W4_OE_vs_Ctrl"
 [7] "condition_O2_KO_vs_Ctrl"   "condition_O2W4_KO_vs_Ctrl"
 [9] "condition_W4_KO_vs_Ctrl"   "condition_W4_OE_vs_Ctrl"  

2 Wald Tests - all

In this section, we export all the pairwise comparisons (W4_OE vs Ctrl, O1_OE vs Ctrl, etc) as an excel sheet. Each sheet will contain one comparison.

Code
# Function to extract results from a given comparison 
get_res <- function(x) {
  res = results(wald_dds, name = x) 
  res$ensembl_gene_id = rownames(res)
  res$external_gene_name = ensembl.genes$external_gene_name[match(rownames(res), ensembl.genes$ensembl_gene_id)]
  res$gene_biotype = ensembl.genes$gene_biotype[match(rownames(res), ensembl.genes$ensembl_gene_id)]
  
  res <- as_tibble(res) %>% 
    dplyr::select(ensembl_gene_id, log2FoldChange, pvalue, padj, gene_biotype, external_gene_name)
  
  return(res)
  }

wald_res_all <- lapply(res_names, get_res)
names(wald_res_all) <- res_names

head(wald_res_all)
$Intercept
# A tibble: 11,145 × 6
   ensembl_gene_id log2FoldChange    pvalue      padj gene_biotype  
   <chr>                    <dbl>     <dbl>     <dbl> <chr>         
 1 FBgn0000003               8.07 2.13e-134 2.71e-134 ncRNA         
 2 FBgn0000008              11.4  0         0         protein_coding
 3 FBgn0000014              12.6  0         0         protein_coding
 4 FBgn0000015              11.9  0         0         protein_coding
 5 FBgn0000017              12.7  0         0         protein_coding
 6 FBgn0000018               7.91 0         0         protein_coding
 7 FBgn0000022               8.17 0         0         protein_coding
 8 FBgn0000024               9.42 0         0         protein_coding
 9 FBgn0000028               7.81 0         0         protein_coding
10 FBgn0000032              10.3  0         0         protein_coding
# ℹ 11,135 more rows
# ℹ 1 more variable: external_gene_name <chr>

$condition_O1_KO_vs_Ctrl
# A tibble: 11,145 × 6
   ensembl_gene_id log2FoldChange    pvalue     padj gene_biotype  
   <chr>                    <dbl>     <dbl>    <dbl> <chr>         
 1 FBgn0000003             -0.323 0.486     0.604    ncRNA         
 2 FBgn0000008             -0.464 0.000288  0.00244  protein_coding
 3 FBgn0000014             -0.548 0.00233   0.0112   protein_coding
 4 FBgn0000015             -0.339 0.0353    0.0810   protein_coding
 5 FBgn0000017              0.154 0.280     0.398    protein_coding
 6 FBgn0000018              0.323 0.0848    0.158    protein_coding
 7 FBgn0000022             -0.124 0.658     0.751    protein_coding
 8 FBgn0000024             -1.17  0.0000568 0.000739 protein_coding
 9 FBgn0000028             -0.185 0.515     0.629    protein_coding
10 FBgn0000032              0.858 0.000360  0.00291  protein_coding
# ℹ 11,135 more rows
# ℹ 1 more variable: external_gene_name <chr>

$condition_O1_OE_vs_Ctrl
# A tibble: 11,145 × 6
   ensembl_gene_id log2FoldChange  pvalue   padj gene_biotype external_gene_name
   <chr>                    <dbl>   <dbl>  <dbl> <chr>        <chr>             
 1 FBgn0000003             -0.250 0.590   0.724  ncRNA        7SLRNA:CR32864    
 2 FBgn0000008             -0.147 0.250   0.429  protein_cod… a                 
 3 FBgn0000014             -0.225 0.211   0.390  protein_cod… abd-A             
 4 FBgn0000015             -0.245 0.127   0.295  protein_cod… Abd-B             
 5 FBgn0000017              0.184 0.194   0.372  protein_cod… Abl               
 6 FBgn0000018              0.208 0.268   0.446  protein_cod… abo               
 7 FBgn0000022              0.217 0.435   0.601  protein_cod… ac                
 8 FBgn0000024             -0.532 0.0653  0.207  protein_cod… Ace               
 9 FBgn0000028             -0.440 0.123   0.288  protein_cod… acj6              
10 FBgn0000032              0.666 0.00565 0.0694 protein_cod… Acph-1            
# ℹ 11,135 more rows

$condition_O1O2_KO_vs_Ctrl
# A tibble: 11,145 × 6
   ensembl_gene_id log2FoldChange      pvalue       padj gene_biotype  
   <chr>                    <dbl>       <dbl>      <dbl> <chr>         
 1 FBgn0000003            0.173   0.709       0.788      ncRNA         
 2 FBgn0000008           -0.343   0.00724     0.0229     protein_coding
 3 FBgn0000014           -0.678   0.000167    0.00123    protein_coding
 4 FBgn0000015           -0.477   0.00309     0.0119     protein_coding
 5 FBgn0000017            0.214   0.131       0.219      protein_coding
 6 FBgn0000018            0.00347 0.985       0.991      protein_coding
 7 FBgn0000022           -0.551   0.0500      0.104      protein_coding
 8 FBgn0000024           -1.49    0.000000336 0.00000930 protein_coding
 9 FBgn0000028           -0.309   0.276       0.391      protein_coding
10 FBgn0000032            0.893   0.000203    0.00143    protein_coding
# ℹ 11,135 more rows
# ℹ 1 more variable: external_gene_name <chr>

$condition_O1W4_KO_vs_Ctrl
# A tibble: 11,145 × 6
   ensembl_gene_id log2FoldChange   pvalue     padj gene_biotype  
   <chr>                    <dbl>    <dbl>    <dbl> <chr>         
 1 FBgn0000003              0.510 2.70e- 1 3.47e- 1 ncRNA         
 2 FBgn0000008             -0.889 5.53e-12 1.26e-10 protein_coding
 3 FBgn0000014             -1.11  7.27e-10 9.37e- 9 protein_coding
 4 FBgn0000015             -0.876 6.17e- 8 4.62e- 7 protein_coding
 5 FBgn0000017              0.756 1.02e- 7 7.23e- 7 protein_coding
 6 FBgn0000018              0.238 2.09e- 1 2.80e- 1 protein_coding
 7 FBgn0000022             -0.447 1.13e- 1 1.67e- 1 protein_coding
 8 FBgn0000024             -1.30  8.41e- 6 3.76e- 5 protein_coding
 9 FBgn0000028             -0.813 5.22e- 3 1.16e- 2 protein_coding
10 FBgn0000032              1.77  1.46e-13 5.19e-12 protein_coding
# ℹ 11,135 more rows
# ℹ 1 more variable: external_gene_name <chr>

$condition_O1W4_OE_vs_Ctrl
# A tibble: 11,145 × 6
   ensembl_gene_id log2FoldChange  pvalue   padj gene_biotype external_gene_name
   <chr>                    <dbl>   <dbl>  <dbl> <chr>        <chr>             
 1 FBgn0000003             0.269  0.561   0.773  ncRNA        7SLRNA:CR32864    
 2 FBgn0000008            -0.213  0.0944  0.301  protein_cod… a                 
 3 FBgn0000014            -0.314  0.0805  0.279  protein_cod… abd-A             
 4 FBgn0000015            -0.323  0.0450  0.209  protein_cod… Abd-B             
 5 FBgn0000017            -0.116  0.413   0.663  protein_cod… Abl               
 6 FBgn0000018             0.117  0.532   0.754  protein_cod… abo               
 7 FBgn0000022             0.473  0.0865  0.288  protein_cod… ac                
 8 FBgn0000024            -0.900  0.00185 0.0364 protein_cod… Ace               
 9 FBgn0000028             0.0268 0.924   0.971  protein_cod… acj6              
10 FBgn0000032             0.136  0.571   0.779  protein_cod… Acph-1            
# ℹ 11,135 more rows
Code
# export this
openxlsx::write.xlsx(wald_res_all, 
                     file = paste0(output_dir, "logfoldchanges.xlsx"), 
                     sheetName = names(wald_res_all), 
                     rowNames = FALSE)

3 Wald Tests

3.1 W4_OE

Code
### Note the order of W4_OE and CTRL - W4_OE here is the numerator, CTRL is the denominator

results_W4OE_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "W4_OE", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_W4_OE_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_W4OE_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 47
Code
### MA plot
generate_ma(res)

Code
### Volcano plot
generate_volcano(res, "W4_OE", "Ctrl")

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 26
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "7 enrichments found"
Code
ego
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0046058 cAMP metabolic process 2/19 11/8880 0.0002358 0.0337452 0.0274482 CG32301/CG32305 2
GO:0050878 regulation of body fluid levels 2/19 14/8880 0.0003887 0.0337452 0.0274482 Spn28Da/CG42259 2
GO:0007193 adenylate cyclase-inhibiting G protein-coupled receptor signaling pathway 2/19 16/8880 0.0005113 0.0337452 0.0274482 CG32301/CG32305 2
GO:0052652 cyclic purine nucleotide metabolic process 2/19 20/8880 0.0008054 0.0352071 0.0286374 CG32301/CG32305 2
GO:0009190 cyclic nucleotide biosynthetic process 2/19 21/8880 0.0008891 0.0352071 0.0286374 CG32301/CG32305 2
GO:0009187 cyclic nucleotide metabolic process 2/19 26/8880 0.0013672 0.0451174 0.0366983 CG32301/CG32305 2
GO:0007189 adenylate cyclase-activating G protein-coupled receptor signaling pathway 2/19 37/8880 0.0027627 0.0781446 0.0635626 CG32301/CG32305 2
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 83
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "4 enrichments found"
Code
ego
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0044242 cellular lipid catabolic process 5/55 86/8880 0.0001802 0.0654606 0.0596368 Pepck2/CG7900/CG15533/CG15534/CG3699 5
GO:0046395 carboxylic acid catabolic process 5/55 106/8880 0.0004774 0.0654606 0.0596368 cn/Oat/Pepck2/CG7900/CG3699 5
GO:0016054 organic acid catabolic process 5/55 107/8880 0.0004984 0.0654606 0.0596368 cn/Oat/Pepck2/CG7900/CG3699 5
GO:0016042 lipid catabolic process 5/55 117/8880 0.0007496 0.0738387 0.0672695 Pepck2/CG7900/CG15533/CG15534/CG3699 5
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster
Reading KEGG annotation online: "https://rest.kegg.jp/link/dme/pathway"...
Reading KEGG annotation online: "https://rest.kegg.jp/list/pathway/dme"...
Reading KEGG annotation online: "https://rest.kegg.jp/conv/ncbi-geneid/dme"...

[1] "2 enrichments found"
Code
write.csv(up_ekegg, file = paste0(output_dir, "w4oe_vs_ctrl_up_kegg.csv"))
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "10 enrichments found"
Code
write.csv(down_ekegg, file = paste0(output_dir, "w4oe_vs_ctrl_down_kegg.csv"))

3.2 O1_OE

Code
### Note the order of W4_OE and CTRL - W4_OE here is the numerator, CTRL is the denominator

results_O1OE_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "O1_OE", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_O1_OE_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_O1OE_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 571
Code
### MA plot
generate_ma(res)

Code
### Volcano plot
generate_volcano(res, "O1_OE", "Ctrl")
Warning: Removed 217 rows containing missing values (`geom_point()`).

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 325
Code
ego_up <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "168 enrichments found"
Code
#ego_up
Code
interesting_pathways <- c("nuclear division", 
                          "regulation of cell cycle",
                          "blastoderm segmentation",
                          "embryonic pattern specification",
                          "anterior/posterior axis specification", 
                          "nuclear chromosome segregation", 
                          "mitotic cell cycle",
                          "DNA repair", 
                          "mismatch repair",
                          "spindle organization")

# Takes in an ego tibble, subsets interesting pathways. 
custom_ego_table(ego_up, interesting_pathways)
Description GeneRatio p.adjust geneID
nuclear division 33/275 0.0000000 bru1/aub/Axs/cort/cup/CycA/CycB/Grip91/mus301/spn-E/tws/yem/gammaTub37C/dhd/mst/Grip128/Ercc1/c(2)M/APC7/Nek2/Fign/Grip71/mms4/Mos/Sse/CG6985/mus81/cid/Nnf1b/Zwilch/alphaTub67C/gwl/vas
anterior/posterior axis specification 21/275 0.0000002 aub/bcd/cad/cup/esc/exu/mus301/nos/spn-E/swa/tor/orb/fs(1)M3/Grip128/krimp/CG9925/armi/clos/vret/qin/vas
nuclear chromosome segregation 23/275 0.0000009 aub/Axs/cort/cup/CycA/CycB/mus301/spn-E/tws/gammaTub37C/mst/Ercc1/c(2)M/APC7/Fign/mms4/Mos/Sse/mus81/cid/Nnf1b/Zwilch/vas
mitotic cell cycle 39/275 0.0000012 bru1/aub/CycA/CycB/png/fs(1)Ya/gnu/Grip91/mus101/spn-E/tws/gammaTub37C/CycJ/l(2)dtl/lok/mre11/mst/Rad9/Grip128/ana2/APC7/Nek2/Fign/Grip71/mms4/Ctf4/nopo/Sse/mus81/cid/Nnf1b/Tap42/Zwilch/del/alphaTub67C/gwl/grp/Cdk7/vas
DNA repair 23/275 0.0000014 Pms2/Irbp/spel1/tos/lok/mre11/agt/Rad9/Ercc1/Xrcc2/CG13690/Nse1/mms4/Ctf4/rad50/Msh6/Snm1/CG10445/CG9601/timeout/CG6985/mus81/gkt
mismatch repair 6/275 0.0000543 Pms2/spel1/tos/Ercc1/CG13690/Msh6
regulation of cell cycle 25/275 0.0001969 bru1/CycA/CycB/png/gnu/mus101/tws/CycJ/l(2)dtl/lok/mre11/Rad9/APC7/Nek2/Grip71/mms4/Mos/nopo/timeout/mus81/armi/Zwilch/gwl/grp/Cdk7
blastoderm segmentation 19/275 0.0002392 aub/bcd/cad/cup/exu/nos/spn-E/swa/tor/orb/fs(1)M3/Grip128/krimp/CG9925/armi/clos/vret/qin/vas
embryonic pattern specification 19/275 0.0006195 aub/bcd/cad/cup/exu/nos/spn-E/swa/tor/orb/fs(1)M3/Grip128/krimp/CG9925/armi/clos/vret/qin/vas
spindle organization 13/275 0.0018714 Axs/Grip91/gammaTub37C/mst/Grip128/spd-2/Grip71/CG18190/cid/Tap42/del/alphaTub67C/grp
Code
custom_ego(ego_up, interesting_pathways,
           title = "GO, O1_OE vs Ctrl")

Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 1026
Code
ego_down <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "16 enrichments found"
Code
interesting_pathways <- c("DNA replication", 
                          "DNA repair",
                          "recombinational repair",
                          "double-strand break repair")

custom_ego(ego_down, interesting_pathways, 
           title = "GO, downregulated, O1OE vs Ctrl")

Code
custom_ego_table(ego_down, interesting_pathways)
Description GeneRatio p.adjust geneID
DNA replication 33/819 0.0000068 Blm/otu/pn/spn-A/Cdk1/tam/PCNA/CycE/RnrL/DNApol-alpha50/Orc5/Mcm5/Orc6/Cdc45/Pol31/Cdc7/RfC38/CG4078/Psf3/DNAlig4/Pif1/Marcal1/CG3430/RfC3/PCNA2/Mcm10/CG7457/RecQ4/Chrac-14/RfC4/Psf2/DNAlig1/DNApol-delta
DNA repair 39/819 0.0004417 aurA/lds/DNApol-zeta/mus304/Blm/phr/spn-A/spn-B/tam/PCNA/Parp/Mcm5/CG14814/Rad17/Cdc45/Rad23/Cdc7/RfC38/CG4078/DNAlig4/Pif1/Marcal1/Mnn1/CG5181/RfC3/Nse4/PCNA2/CG8679/CG9062/CG5543/CG7457/dgrn/CG5466/RecQ4/RfC4/Xpd/Psf2/DNAlig1/DNApol-delta
recombinational repair 17/819 0.0115503 DNApol-zeta/Blm/spn-A/spn-B/Mcm5/Cdc45/Cdc7/CG4078/Pif1/Marcal1/CG5181/CG8679/CG9062/CG7457/dgrn/RecQ4/Psf2
double-strand break repair 20/819 0.0938205 DNApol-zeta/Blm/spn-A/spn-B/Parp/Mcm5/Cdc45/Cdc7/CG4078/DNAlig4/Pif1/Marcal1/CG5181/CG8679/CG9062/CG5543/CG7457/dgrn/RecQ4/Psf2
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "6 enrichments found"
Code
custom_ekegg(up_ekegg, interesting_pathways = up_ekegg$Description,
             title = "KEGG, up, O1OE vs Ctrl")

Code
write.csv(up_ekegg, file = paste0(output_dir, "o1oe_vs_ctrl_up_kegg.csv"))
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "16 enrichments found"
Code
down_ekegg
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme03030 dme03030 DNA replication - Drosophila melanogaster (fruit fly) 11/273 34/3488 0.0000310 0.0027965 0.0023493 PCNA/DNApol-alpha50/Mcm5/Pol31/RfC38/RfC3/PCNA2/Chrac-14/CG30105/RfC4/DNApol-delta 11
dme03410 dme03410 Base excision repair - Drosophila melanogaster (fruit fly) 10/273 30/3488 0.0000538 0.0027965 0.0023493 tam/PCNA/Parp/Pol31/RfC38/RfC3/PCNA2/Chrac-14/RfC4/DNApol-delta 10
dme00983 dme00983 Drug metabolism - other enzymes - Drosophila melanogaster (fruit fly) 19/273 94/3488 0.0000842 0.0029189 0.0024521 bur/GstD7/RnrL/CRMP/Ugt37C1/Ugt36E1/Ugt49B1/nmdyn-D7/GstT3/CG15117/pyd3/Ugt303A1/Ugt37D1/Ugt37C2/GstD10/GstE5/GstE2/su(r)/dUTPase 19
dme03430 dme03430 Mismatch repair - Drosophila melanogaster (fruit fly) 7/273 19/3488 0.0003695 0.0096076 0.0080711 PCNA/Pol31/RfC38/RfC3/PCNA2/RfC4/DNApol-delta 7
dme01240 dme01240 Biosynthesis of cofactors - Drosophila melanogaster (fruit fly) 21/273 138/3488 0.0020569 0.0427837 0.0359418 Dhod/Alp4/Ugt37C1/Ugt36E1/Ugt49B1/nmdyn-D7/Nadsyn/CG12237/Naprt/CG15117/Alp7/Alp8/Alp10/Alp9/UGP/Hpd/Ugt303A1/Ugt37D1/Ugt37C2/sgll/Alp1 21
dme00770 dme00770 Pantothenate and CoA biosynthesis - Drosophila melanogaster (fruit fly) 5/273 14/3488 0.0031444 0.0481389 0.0404406 b/CRMP/CG5618/pyd3/su(r) 5
dme03013 dme03013 Nucleocytoplasmic transport - Drosophila melanogaster (fruit fly) 16/273 99/3488 0.0037598 0.0481389 0.0404406 RanGAP/msk/Nup93-1/Ranbp21/Tnpo-SR/Nup54/Rae1/thoc7/thoc6/Nup43/Nup58/Nup133/Ndc1/Arts/Fs(2)Ket/Kr-h2 16
dme00730 dme00730 Thiamine metabolism - Drosophila melanogaster (fruit fly) 6/273 21/3488 0.0043171 0.0481389 0.0404406 Alp4/Alp7/Alp8/Alp10/Alp9/Alp1 6
dme00982 dme00982 Drug metabolism - cytochrome P450 - Drosophila melanogaster (fruit fly) 12/273 66/3488 0.0044156 0.0481389 0.0404406 GstD7/Fdh/Ugt37C1/Ugt36E1/Ugt49B1/GstT3/Ugt303A1/Ugt37D1/Ugt37C2/GstD10/GstE5/GstE2 12
dme00053 dme00053 Ascorbate and aldarate metabolism - Drosophila melanogaster (fruit fly) 8/273 35/3488 0.0046287 0.0481389 0.0404406 Ugt37C1/Ugt36E1/Ugt49B1/CG15117/CG6910/Ugt303A1/Ugt37D1/Ugt37C2 8
dme03420 dme03420 Nucleotide excision repair - Drosophila melanogaster (fruit fly) 10/273 51/3488 0.0052402 0.0491647 0.0413024 PCNA/Rad23/Pol31/RfC38/RfC3/PCNA2/Chrac-14/RfC4/Xpd/DNApol-delta 10
dme00980 dme00980 Metabolism of xenobiotics by cytochrome P450 - Drosophila melanogaster (fruit fly) 12/273 68/3488 0.0056729 0.0491647 0.0413024 GstD7/Fdh/Ugt37C1/Ugt36E1/Ugt49B1/GstT3/Ugt303A1/Ugt37D1/Ugt37C2/GstD10/GstE5/GstE2 12
dme00240 dme00240 Pyrimidine metabolism - Drosophila melanogaster (fruit fly) 8/273 38/3488 0.0078333 0.0626666 0.0526450 Dhod/RnrL/dnk/CRMP/nmdyn-D7/pyd3/su(r)/dUTPase 8
dme00410 dme00410 beta-Alanine metabolism - Drosophila melanogaster (fruit fly) 6/273 25/3488 0.0108042 0.0802600 0.0674249 b/CRMP/CG17896/CG5618/pyd3/su(r) 6
dme00830 dme00830 Retinol metabolism - Drosophila melanogaster (fruit fly) 7/273 33/3488 0.0121210 0.0840390 0.0705996 Fdh/Ugt37C1/Ugt36E1/Ugt49B1/Ugt303A1/Ugt37D1/Ugt37C2 7
dme00040 dme00040 Pentose and glucuronate interconversions - Drosophila melanogaster (fruit fly) 9/273 50/3488 0.0139376 0.0905942 0.0761065 Ugt37C1/Ugt36E1/Ugt49B1/Had2/CG15117/UGP/Ugt303A1/Ugt37D1/Ugt37C2 9
Code
interesting_pathways = c("DNA replication",
                         "Base excision repair",
                         "Mismatch repair",
                         "Nucleotide excision repair")

custom_ekegg(down_ekegg, 
             interesting_pathways, 
             title = "KEGG, down, O1OE vs Ctrl")

Code
write.csv(down_ekegg, file = paste0(output_dir, "o1oe_vs_ctrl_down_kegg.csv"))

3.3 O1W4_OE

Code
results_O1W4OE_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "O1W4_OE", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_O1W4_OE_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_O1W4OE_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 662
Code
### MA plot
generate_ma(res)
Warning: Removed 1 rows containing missing values (`geom_point()`).

Code
### Volcano plot
generate_volcano(res, "O1W4_OE", "Ctrl")
Warning: Removed 1 rows containing missing values (`geom_point()`).

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 141
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "15 enrichments found"
Code
nrow(ego)
[1] 15
Code
write.csv(ego, file = paste0(output_dir, "o1w4oe_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0098813 nuclear chromosome segregation 10/103 178/8880 3.83e-05 0.0128926 0.0118573 cort/mus301/gammaTub37C/Ercc1/Fign/RfC3/mms4/Mos/Nnf1a/vas 10
GO:0061982 meiosis I cell cycle process 8/103 110/8880 3.85e-05 0.0128926 0.0118573 cort/mus301/gammaTub37C/Ercc1/Grip71/mms4/del/vas 8
GO:0006281 DNA repair 10/103 186/8880 5.57e-05 0.0128926 0.0118573 aurA/Pms2/tos/Rad9/Ercc1/Xrcc2/Nse1/RfC3/mms4/gkt 10
GO:0048285 organelle fission 13/103 312/8880 6.01e-05 0.0128926 0.0118573 aurA/cort/mus301/gammaTub37C/dhd/Ercc1/Fign/Grip71/mms4/Mos/Nnf1a/CG33474/vas 13
GO:0006974 cellular response to DNA damage stimulus 11/103 233/8880 7.63e-05 0.0128926 0.0118573 aurA/Pms2/tos/Rad9/Ercc1/Xrcc2/Nse1/RfC3/mms4/armi/gkt 11
GO:0007059 chromosome segregation 10/103 198/8880 9.42e-05 0.0128926 0.0118573 cort/mus301/gammaTub37C/Ercc1/Fign/RfC3/mms4/Mos/Nnf1a/vas 10
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 1063
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "42 enrichments found"
Code
nrow(ego)
[1] 42
Code
write.csv(ego, file = paste0(output_dir, "o1w4oe_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0040003 chitin-based cuticle development 37/809 147/8880 0.0e+00 0.0000148 0.0000145 Gld/Ccp84Ag/fz2/Gasp/Idgf4/Vajk3/Vajk1/obst-A/TwdlE/Vajk2/Twdlbeta/Cpr49Ae/Cpr49Ag/Cpr50Cb/Cpr51A/Cpr62Bb/Cht7/Cpr64Ab/Cpr64Ad/Cpr67B/Cpr67Fa1/Cpr78Cc/TwdlW/TwdlO/TwdlD/TwdlQ/Vajk4/TwdlH/fred/Cpr66D/TwdlY/TwdlX/dpy/Tb/pot/drd/verm 37
GO:0042335 cuticle development 40/809 185/8880 2.0e-07 0.0001548 0.0001517 Gld/Ccp84Ag/Cyp4g1/fz2/Gasp/Idgf4/Vajk3/Vajk1/obst-A/snsl/TwdlE/Vajk2/Twdlbeta/Cpr49Ae/Cpr49Ag/Cpr50Cb/Cpr51A/Cpr62Bb/Cht7/Cpr64Ab/Cpr64Ad/Cpr67B/Cpr67Fa1/Cpr78Cc/TwdlW/TwdlO/TwdlD/TwdlQ/snu/Vajk4/TwdlH/fred/Cpr66D/TwdlY/TwdlX/dpy/Tb/pot/drd/verm 40
GO:0003008 system process 77/809 472/8880 2.0e-07 0.0001548 0.0001517 nAChRalpha1/nAChRbeta1/Col4a1/eag/Hn/Mlc2/na/Pkc98E/Rh4/Sh/ss/w/Tm2/5-HT1A/up/ct/Rdl/me31B/Shal/Galphaf/Mdr50/Dop1R1/Drip/H15/nompC/mGluR/ppk/b6/GluClalpha/Scgdelta/CG4928/nAChRalpha5/Nha1/Obp44a/KCNQ/CAP/dpr13/Eglp2/DopEcR/Or71a/CG14636/Nrx-1/Root/robo3/Obp83g/ClC-a/wry/dpr8/dpr11/Nmdar2/Dop2R/Unc-89/dpr12/Nox/sif/futsch/brp/drd/Mp/2mit/stj/rdgA/Npc1b/Shab/cac/Fic/dysc/CG43778/Ca-alpha1T/Mhc/tn/rad/rdgC/Fhos/SclB/Tmc/para 77
GO:0003012 muscle system process 13/809 33/8880 2.8e-06 0.0015546 0.0015231 Mlc2/Tm2/up/CG4928/KCNQ/robo3/Unc-89/Nox/Mp/Mhc/tn/Fhos/SclB 13
GO:0099536 synaptic signaling 51/809 296/8880 5.3e-06 0.0022051 0.0021605 Ace/nAChRalpha1/mAChR-A/nAChRbeta1/Cat/na/Sh/slo/5-HT1A/Rdl/Cbp53E/Shal/Rab3/Dop1R1/nSyb/stnB/mGluR/b6/GluClalpha/unc-13/GABA-B-R2/nAChRalpha5/Rab3-GEF/stas/Nlg2/CG6583/atilla/Dyb/spz5/unc-13-4A/twf/Nrx-1/Tusp/CDase/dj-1beta/kek6/CG31676/scramb1/Rim/Nlg3/brp/Dys/stj/Shab/cac/Ih/nwk/cpo/dysc/Fife/Esyt2 51
GO:0007268 chemical synaptic transmission 50/809 291/8880 7.1e-06 0.0022051 0.0021605 Ace/nAChRalpha1/mAChR-A/nAChRbeta1/Cat/na/Sh/slo/5-HT1A/Rdl/Cbp53E/Shal/Rab3/Dop1R1/nSyb/stnB/mGluR/b6/GluClalpha/unc-13/GABA-B-R2/nAChRalpha5/Rab3-GEF/stas/Nlg2/CG6583/atilla/Dyb/spz5/unc-13-4A/twf/Nrx-1/CDase/dj-1beta/kek6/CG31676/scramb1/Rim/Nlg3/brp/Dys/stj/Shab/cac/Ih/nwk/cpo/dysc/Fife/Esyt2 50
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "2 enrichments found"
Code
nrow(up_ekegg)
[1] 2
Code
write.csv(up_ekegg, file = paste0(output_dir, "o1w4oe_vs_ctrl_up_kegg.csv"))
Code
head(up_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme03430 dme03430 Mismatch repair - Drosophila melanogaster (fruit fly) 3/37 19/3488 0.0009476 0.0397980 0.0309208 Pms2/tos/RfC3 3
dme03460 dme03460 Fanconi anemia pathway - Drosophila melanogaster (fruit fly) 3/37 26/3488 0.0024158 0.0507308 0.0394149 Pms2/Ercc1/mms4 3
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "9 enrichments found"
Code
nrow(down_ekegg)
[1] 9
Code
write.csv(down_ekegg, file = paste0(output_dir, "o1w4oe_vs_ctrl_down_kegg.csv"))
Code
head(down_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme04080 dme04080 Neuroactive ligand-receptor interaction - Drosophila melanogaster (fruit fly) 17/206 67/3488 0.0000002 0.0000165 0.0000143 mAChR-A/5-HT1A/Rdl/deltaTry/gammaTry/epsilonTry/Dop1R1/mGluR/GABA-B-R2/CCKLR-17D3/GABA-B-R3/SIFaR/lambdaTry/kappaTry/Nmdar2/Dop2R/GABA-B-R1 17
dme00790 dme00790 Folate biosynthesis - Drosophila melanogaster (fruit fly) 10/206 40/3488 0.0000737 0.0037221 0.0032198 Hn/Pcd/CG12116/Alp6/Alp7/Alp8/Alp10/GstO3/Alp5/Mocs1 10
dme04512 dme04512 ECM-receptor interaction - Drosophila melanogaster (fruit fly) 5/206 11/3488 0.0002361 0.0079484 0.0068756 Col4a1/vkg/Hml/CG3168/Tsp 5
dme00230 dme00230 Purine metabolism - Drosophila melanogaster (fruit fly) 13/206 92/3488 0.0025482 0.0455990 0.0394447 ry/Ac76E/Gycalpha99B/dnk/nmdyn-D7/nmdyn-D6/Pde6/Gyc88E/CG34357/DhpD/Pde1c/Pde8/AOX1 13
dme01240 dme01240 Biosynthesis of cofactors - Drosophila melanogaster (fruit fly) 17/206 138/3488 0.0027009 0.0455990 0.0394447 v/Ugt37B1/nmdyn-D7/CG12116/nmdyn-D6/Alp6/Alp7/Alp8/Alp10/Hpd/Alp5/Ugt303A1/Ugt35E2/Ugt37D1/Ugt37C2/Mocs1/CG44243 17
dme00410 dme00410 beta-Alanine metabolism - Drosophila melanogaster (fruit fly) 6/206 25/3488 0.0027089 0.0455990 0.0394447 b/Aldh-III/CG9527/CG15093/Gabat/su(r) 6

3.4 O2_KO

Code
results_O2KO_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "O2_KO", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_O2_KO_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_O2KO_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 320
Code
### MA plot
generate_ma(res)

Code
### Volcano plot
generate_volcano(res, "O2_KO", "Ctrl")

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 129
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster
Warning in plotEGO_dm(up_genes, universe = rownames(dds), ont = "BP", title =
"GO, upregulated genes"): No GO enrichment found. Returning a NA tibble.
Code
nrow(ego)
[1] 1
Code
write.csv(ego, file = paste0(output_dir, "O2KO_vs_ctrl_down_go.csv"))
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 500
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster
Warning in plotEGO_dm(down_genes, universe = rownames(dds), ont = "BP", : No GO
enrichment found. Returning a NA tibble.
Code
nrow(ego)
[1] 1
Code
write.csv(ego, file = paste0(output_dir, "O2KO_vs_ctrl_down_go.csv"))
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "2 enrichments found"
Code
nrow(up_ekegg)
[1] 2
Code
write.csv(up_ekegg, file = paste0(output_dir, "O2KO_vs_ctrl_up_kegg.csv"))
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "3 enrichments found"
Code
nrow(down_ekegg)
[1] 3
Code
write.csv(down_ekegg, file = paste0(output_dir, "O2KO_vs_ctrl_down_kegg.csv"))

3.5 W4_KO

Code
results_W4KO_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "W4_KO", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_W4_KO_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_W4KO_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 4355
Code
### MA plot
generate_ma(res)

Code
### Volcano plot
generate_volcano(res, "W4_KO", "Ctrl")

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 2149
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "398 enrichments found"
Code
nrow(ego)
[1] 398
Code
write.csv(ego, file = paste0(output_dir, "W4KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0051276 chromosome organization 141/1819 298/8880 0 0 0 asp/aub/Cam/CycA/CycB/Hmr/kel/pho/mei-S332/mod(mdg4)/mr/ms(3)K81/mu2/mus301/mus304/Blm/ncd/pim/polo/spn-B/spn-D/spn-E/sub/Su(var)3-7/Su(var)205/Su(var)2-10/thr/Klp61F/Klp67A/shtd/rhi/wapl/piwi/enc/tws/Top1/gammaTub37C/Parp/pav/Irbp/Mtor/prod/glu/Cul1/spel1/HDAC1/pgc/Stat92E/mre11/Mcm7/chb/Nup154/DCTN2-p50/san/aurB/SuUR/Bub3/Bap60/Cdc16/Mcm6/Rok/cav/Mad1/Nipped-B/Su(var)2-HP2/RecQ5/MRG15/Nup93-1/SMC2/ball/MBD-like/RfC38/CG4078/APC7/Fign/Pif1/HP6/Marcal1/lid/Nuf2/borr/Mcm10/Nup50/Mos/mars/tej/krimp/Nup62/Oxp/Panx/Rae1/Sse/mad2/nudE/PCID2/Smn/Snm1/Hat1/CG8223/BigH1/Mau2/cal1/Nup133/HP1c/Nup98-96/Mink/Kat60/cana/Ku80/Chrac-14/SA-2/Chro/tefu/Haspin/Rif1/nenya/Nnf1b/His1:CG33801/His1:CG33861/Zwilch/del/tum/egg/Spc25/lmgA/pds5/Su(var)3-3/gwl/gammaTub23C/vtd/Incenp/zuc/CoRest/Xpd/Psf2/Mi-2/ver/BubR1/Topors/solo/Top2 141
GO:0006259 DNA metabolic process 141/1819 329/8880 0 0 0 ben/chif/DNaseII/dup/gnu/hay/Hmr/pho/mago/ms(3)K81/mu2/mus101/DNApol-zeta/mus304/Blm/mus312/otu/plu/pn/stwl/spn-B/spn-D/spn-E/Su(var)205/Su(var)2-10/Cdk1/tam/Xpc/RpII15/piwi/Top1/PCNA/Parp/Pms2/RnrL/Irbp/mof/CkIalpha/spel1/tos/smid/mre11/Mcm7/rnh1/Pp4-19C/Orc6/Parg/CG14803/CG14805/flfl/Drep1/agt/SuUR/Rad9/Rad17/Mcm6/cav/Nipped-B/Dek/ADD1/CG6171/XRCC1/Rad23/RecQ5/Pol31/Ercc1/c(2)M/RfC38/CG4078/Caf1-180/cerv/Arp8/CG13690/Tfb4/Pif1/HP6/Marcal1/Nse1/Rat1/Mnn1/CG5181/Nse4/CG10336/CG9272/Mcm10/CG8679/geminin/CG9062/Dcr-2/mip40/DMAP1/pita/Psf1/Uev1A/Sse/Cdc6/CG3448/ssp/Msh6/Ssl1/Snm1/dgrn/MAGE/CG9601/RnpS1/Tctp/timeout/BigH1/CG5466/CG6985/Sld5/mus81/ago/Ku80/Chrac-14/tefu/Rif1/nenya/pch2/Tengl3/CG31807/His1:CG33801/His1:CG33861/AlkB/Ino80/CG42374/pds5/gkt/pic/mei-38/vtd/Prp19/Xpd/Psf2/ver/Nup160/Gen/Caf1-55/Tfb5/RPA3/Top2 141
GO:0000278 mitotic cell cycle 178/1819 465/8880 0 0 0 bru1/asp/aub/awd/ben/csw/CycA/CycB/dsh/png/dup/fzy/gnu/Hmr/Grip91/l(3)mbt/sti/Rcc1/Map205/twe/mor/mr/mud/mus101/mus304/ncd/pim/plu/polo/PpV/pum/smo/spn-E/sub/Su(var)205/thr/alphaTub84D/Cdk1/Klp61F/Klp67A/shtd/wapl/Su(H)/RpII15/tws/PCNA/PpD3/gammaTub37C/CycD/CycJ/dia/Mad/pav/tsr/Wee1/Arp53D/l(2)dtl/pnut/Mtor/prod/noi/glu/CycB3/Ran/mre11/asun/chb/DCTN2-p50/Pp4-19C/Orc6/san/aurB/Dlc90F/Bub3/Cdc16/Rad9/Rad17/Mcm6/JTBR/Rok/Mad1/Nipped-B/Su(var)2-HP2/Grip75/Grip163/Apc2/RecQ5/SMC2/ball/Pol31/CHORD/Aos1/APC7/Nek2/Aladin/PPP4R2r/Sap30/CG8173/Rab35/Fign/Mnn1/Nuf2/CG5181/borr/CG10336/Grip71/Mcm10/geminin/dgt5/mars/dgt3/Rae1/JMJD5/Psf1/Sse/mad2/Cdc6/nudE/Syx13/Skp2/PEK/Tctp/Rrp6/Mau2/cal1/CG13599/Dis3/Mink/Sld5/CG5003/dgt6/pont/pins/Kat60/cana/MED22/mus81/ago/Chro/tefu/Haspin/Rif1/Nnf1b/Tap42/Claspin/Zwilch/del/ALiX/tum/Spc25/alphaTub67C/lmgA/zld/bora/pds5/gwl/gammaTub23C/wisp/mei-38/vtd/Incenp/grp/scra/nej/SmD2/Xpd/aPKC/CG43295/Cdk7/BubR1/Caf1-55/Det/smt3/Smr/Ptp61F/Tlk/Top2/betaTub56D 178
GO:0051726 regulation of cell cycle 127/1819 303/8880 0 0 0 bru1/ben/chif/CycA/CycB/png/dup/gnu/l(3)mbt/sti/Rcc1/twe/mor/mr/mus101/mus304/ncd/otu/plu/polo/pum/smo/sub/Cdk1/shtd/wapl/Su(H)/tws/CycD/CycJ/mtrm/Mad/Wee1/l(2)dtl/Mtor/mof/Cul1/CycB3/Gcn2/Ran/mre11/asun/chb/Pp4-19C/aurB/SuUR/CG12179/Bub3/Cdc16/Rad9/Rad17/Rok/Tsc1/Mad1/Sep5/CSN5/Stam/RecQ5/CHORD/Aos1/APC7/Nek2/Sap30/CG8173/Hrs/Mnn1/CG5181/CycY/Ing5/CG10336/Grip71/Mcm10/geminin/Br140/dgt5/Mos/mars/dgt3/Rae1/Eaf6/mad2/Cdc6/ssp/Skp2/PEK/CG17816/Tctp/timeout/Hexim/CG5466/CG13599/CG5003/dgt6/cana/mus81/armi/ago/Chro/tefu/Haspin/Rif1/nenya/pch2/Claspin/Zwilch/tum/lmgA/mnb/bora/pds5/gwl/gammaTub23C/Incenp/asl/grp/scra/nej/Xpd/Mi-2/CG43295/Cdk7/BubR1/Caf1-55/SNF4Agamma/Smr/elgi/Tlk 127
GO:1903047 mitotic cell cycle process 133/1819 327/8880 0 0 0 bru1/asp/aub/CycA/CycB/dsh/Hmr/sti/twe/mor/mr/mud/mus101/mus304/ncd/pim/polo/PpV/smo/spn-E/sub/thr/Cdk1/Klp61F/Klp67A/shtd/Su(H)/tws/PCNA/gammaTub37C/CycD/CycJ/dia/Mad/pav/tsr/Wee1/Arp53D/l(2)dtl/pnut/Mtor/prod/glu/CycB3/Ran/mre11/chb/DCTN2-p50/san/aurB/Bub3/Cdc16/Rad9/Rad17/Mcm6/JTBR/Rok/Mad1/Nipped-B/Su(var)2-HP2/Apc2/RecQ5/SMC2/ball/Pol31/CHORD/APC7/Nek2/Aladin/Sap30/CG8173/Rab35/Fign/Mnn1/Nuf2/CG5181/borr/CG10336/dgt5/mars/dgt3/Rae1/JMJD5/Psf1/Sse/mad2/Cdc6/nudE/PEK/Tctp/Mau2/CG13599/Mink/CG5003/dgt6/pont/pins/Kat60/cana/mus81/ago/Chro/tefu/Haspin/Rif1/Nnf1b/Tap42/Claspin/Zwilch/del/ALiX/tum/Spc25/alphaTub67C/lmgA/bora/pds5/gwl/gammaTub23C/wisp/mei-38/vtd/Incenp/grp/scra/nej/Xpd/aPKC/CG43295/BubR1/Det/smt3/Top2 133
GO:0006913 nucleocytoplasmic transport 68/1819 124/8880 0 0 0 Abl/Pka-C1/e(y)2/fu/Rcc1/RanGAP/Su(fu)/Nup214/Ref1/Mtor/Hel25E/Ran/emb/Nup154/Cse1/Cbp80/Dbp80/Kap-alpha1/Tnpo/Torsin/mbo/msk/Kap-alpha3/Nup93-1/Nxt1/Aladin/CG9915/Nup35/CG18259/CG6961/Ranbp21/Ntf-2/tho2/Tnpo-SR/Ripalpha/Nup50/Gle1/Nup54/Nup62/Mapmodulin/Rae1/thoc5/thoc7/PCID2/thoc6/nxf2/Hpr1/Ranbp9/Nup93-2/Nup133/RanBP3/Nup98-96/Nup358/tamo/Arts/Apl/Ranbp16/Nup153/Snup/Mat89Ba/cdm/aPKC/Nup160/Nxf3/Rpb4/ebo/Gp210/Pen 68
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 4433
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "342 enrichments found"
Code
nrow(ego)
[1] 342
Code
write.csv(ego, file = paste0(output_dir, "W4KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0040011 locomotion 251/3620 429/8880 0 0 0 Arl1/brm/bsk/chic/eya/Dab/Lar/E(Pc)/Fas2/Fas3/FER/Glt/grn/if/kay/dlg1/Sos/Gli/Ca-alpha1D/Catsup/Lim3/sas/LanA/robo2/sens/Mlc2/peb/pk/Pu/sm/svp/Egfr/ttk/tup/vir/vn/w/Pp1-87B/ct/Ptp99A/Ptp10D/Ten-m/mew/pros/zfh1/hh/fru/sim/Hrb27C/otk/tok/eIF2beta/ple/robo1/nvy/slbo/Orc3/dac/Con/comm/Tbh/Rac1/Cdc42/htl/Nmdar1/Sdc/trn/Wnt4/tutl/kermit/msn/fbl/Sema1a/Dop1R1/Moe/Tig/nmo/beat-Ia/shot/Cdk5/Gycalpha99B/Rac2/Rho1/Tie/bnl/cher/gcm/mirr/tim/drl/mbc/tap/CadN/Lis-1/NetA/NetB/Rab6/Vang/side/wun/fz2/14-3-3epsilon/ko/dom/Nrk/Dad/CG7179/Pka-R2/caps/Pdf/CG3857/trio/stan/DAAM/Bap55/TBPH/Timp/Cdk5alpha/jbug/bdl/Rich/beat-Ic/beat-Ib/nerfin-1/SoxN/raptor/pod1/Clic/hiw/CG9784/CG6867/wgn/Nup205/daw/Pvr/Bug22/Oli/CarT/Dscam1/Mmp2/ths/Drl-2/ckn/Hr51/unc-5/Syn2/Toll-7/Lrt/spz5/RhoGEF64C/CCDC151/atk/M6/Syn1/laza/ATP6AP2/beat-Va/beat-Vb/Sra-1/beat-IIb/beat-IIa/Lrrk/Ktl/RhoGAP93B/Nrx-1/AP-1sigma/XNP/CG6330/beat-VI/DIP-gamma/RhoGAP100F/Vav/Ephrin/Alk/tai/robo3/comm2/dlp/SCAR/mask/Proc/CG30456/Unc-115a/Trim9/CG32100/Krn/gogo/shep/beat-IIIb/dpr11/Mical/pain/CG34353/Pvf3/spri/lap/stl/vvl/l(3)L1231/iav/sls/psidin/puc/kra/beta-Spec/meigo/mmy/kuz/flr/Diap1/Mp/Pdfr/Dscam3/conv/haf/NT1/wb/pdm3/Bsg/Msp300/mid/Ptpmeg/CadN2/Src64B/Evi5/Acsl/hts/Hmgcr/dally/Stacl/dysc/sli/Sema2b/ab/Mhc/Nrg/CG44153/Neto/beat-IIIa/rdgC/lov/Sec15/Ten-a/unc-104/SelR/ap/br/InR/Sac1/Sema5c/trol 251
GO:0034330 cell junction organization 220/3620 371/8880 0 0 0 Act57B/Appl/arm/baz/Lar/Fas2/Galphas/Gs2/Jra/kay/dlg1/Gli/l(2)gl/LanA/Atpalpha/rut/sax/shi/slo/stau/kst/Ten-m/Fur1/Gad1/Syn/pros/Rap1/gig/Rab3/comm/Rac1/Cdc42/htl/Itgbn/Sdc/cora/sinu/jar/Sema1a/cpb/Nos/repo/Hem/nmo/Cdk5/Nrx-IV/Rab2/Rac2/Rho1/drl/nrv2/Rab11/RtGEF/lbm/wun/fz2/klg/mGluR/Pdk1/Dad/caps/wit/sns/WASp/trio/Tsp2A/Flo1/stan/Syx4/moody/SkpA/Srpk79D/TBPH/wrapper/Top3beta/Tip60/drpr/cindr/kirre/Rich/Vap33/Grip/raptor/pod1/dpr14/Upf2/hiw/Lrp4/dpr18/Frq1/shv/bark/DIP-eta/Nlg2/dpr19/PICK1/Cul2/kune/Tbce/Mmp2/Dyb/Syn2/dpr13/Magi/Liprin-gamma/CG13506/PIP5K59B/CG10339/Mmp1/Eps-15/dpr20/CG13704/crim/Tsf2/CG9628/sff/Dysb/Ssk/cmpy/M6/CG7166/Syn1/beag/ATP6AP2/dpr5/dpr15/Octbeta2R/dpr9/Sra-1/Lrrk/Nrx-1/Octbeta1R/spas/BOD1/DIP-gamma/Cul5/ppk30/RhoGAP100F/mav/Cals/Vav/dbo/Aplip1/dpr1/Prosap/Atx2/dlp/cpx/Fsn/bchs/metro/mesh/Nlg1/dpr17/CG31475/DIP-theta/DIP-zeta/DIP-kappa/dpr10/dpr8/DIP-alpha/Synd/dpr11/Mical/dpr4/dpr3/Frq2/Nlg3/Nlg4/shakB/dpr12/DIP-delta/beta-Spec/futsch/DIP-beta/brp/Mob2/DIP-epsilon/Ca-beta/mtg/p120ctn/Atg1/dpr21/wol/stj/bou/CG42673/Thor/Ank2/CG42748/dpr2/sdt/Rbp/Rab8/Arp3/Src64B/cac/scrib/nwk/dysc/Fife/CaMKII/Nrg/l(3)neo38/CG44325/Neto/CrebB/rg/ringer/Nsf2/GckIII/btsz/Ten-a/unc-104/Mcr/Sac1/wg/Adf1/btl/Imp 220
GO:0003008 system process 269/3620 472/8880 0 0 0 nAChRalpha1/nAChRbeta1/nAChRalpha2/Appl/arm/boss/bsk/bw/salr/Col4a1/eag/eas/EcR/Fas2/Galphas/Hn/if/inaD/Jra/Ca-alpha1D/Mlc1/Mlc2/na/Atpalpha/neur/Pkc98E/ru/rut/Sh/shi/Btk29A/ss/stau/sws/Egfr/vn/w/Pp1-87B/Tm2/5-HT1A/up/ct/Rdl/Ptp10D/Klp64D/Klp68D/mew/Mdr65/Oct-TyrR/Gad1/Syn/pros/gl/TkR99D/drk/svr/hop/14-3-3zeta/Shal/Orc3/Itpr/Galphaf/Lcch3/Tbh/Rac1/Cdc42/Nmdar1/cora/sinu/Obp28a/RyR/Ms/Dop1R1/ine/Mvl/Nos/Fdh/Nrx-IV/cher/tim/drl/ktub/nrv2/Drip/nompA/vri/H15/nompB/Dyrk2/klg/mGluR/ppk/geko/Rfx/CRMP/Clk/btv/teq/CG3078/sns/b6/Oamb/GluClalpha/Scgdelta/Amnionless/Eb1/Inx7/Amph/CG4928/kirre/Kap3/nAChRalpha5/aru/NtR/CG5921/dpr14/Uba5/Tomosyn/hiw/dpr18/Mec2/nAChRbeta3/Ir25a/meng/Nha1/PrBP/nAChRalpha6/dpr19/AstC/MRP/CarT/nrv3/Dscam1/Obp44a/CG2121/KCNQ/CAP/Sans/dpr13/Eglp2/Eglp4/ITP/mri/dpr20/Chd64/DopEcR/Fitm/CG10226/ppk26/Bmcp/Or71a/SecCl/Dysb/Ir76b/SLC22A/CG14636/NPFR/Neurochondrin/ATP6AP2/mura/dpr5/ATP8B/dpr15/dpr9/Nrx-1/CCAP/CG6972/Root/sosie/PQBP1/BOD1/Obp99a/Obp99c/Obp99b/Cad99C/CG15537/G9a/Alk/dpr1/robo3/Socs36E/Atx2/HDAC4/prt/Proc/Obp83g/Obp83ef/Obp83cd/ClC-a/SKIP/dpr17/Gba1b/wry/dpr10/Rbfox1/Krn/dpr8/Cubn/dpr11/pxb/dpr4/Nmdar2/dpr3/Dop2R/Unc-89/pain/Ctr1A/Or19b/Lasp/dbr/dpr12/Nox/sif/Tpi/iav/kra/DCX-EMAP/futsch/brp/Mob2/drd/qvr/Mp/2mit/Atg1/stj/Atg13/pdm3/salm/nemy/Npc1b/Ank2/dpr2/dikar/mid/f/Shab/Tmhs/Src64B/AGO1/Ptr/SERCA/cac/Fic/scrib/Nha2/spoon/dysc/orb2/CG43778/Ca-alpha1T/CaMKII/Mhc/Rgk1/tn/rad/CrebB/rdgC/rg/nAChRalpha4/SclA/SclB/Ggamma30A/Tmc/wupA/Phm/wg/Adf1/para 269
GO:0042330 taxis 200/3620 332/8880 0 0 0 brm/bsk/chic/eya/Dab/Lar/E(Pc)/Fas2/Fas3/FER/Glt/grn/if/dlg1/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/sm/svp/tup/w/Pp1-87B/ct/Ptp99A/Ptp10D/Ten-m/mew/pros/zfh1/hh/fru/otk/tok/eIF2beta/ple/robo1/nvy/dac/Con/comm/Rac1/Cdc42/Nmdar1/Sdc/trn/Wnt4/tutl/kermit/Sema1a/Dop1R1/Moe/Tig/nmo/beat-Ia/shot/Cdk5/Gycalpha99B/Rac2/Rho1/cher/gcm/mirr/tim/drl/tap/CadN/Lis-1/NetA/NetB/Rab6/Vang/side/wun/fz2/14-3-3epsilon/ko/dom/Nrk/Dad/Pka-R2/caps/Pdf/CG3857/trio/stan/DAAM/Bap55/Timp/Cdk5alpha/jbug/bdl/Rich/beat-Ic/beat-Ib/nerfin-1/SoxN/raptor/pod1/Clic/hiw/CG6867/wgn/Nup205/daw/Oli/CarT/Dscam1/Mmp2/Drl-2/ckn/Hr51/unc-5/Toll-7/Lrt/spz5/RhoGEF64C/CCDC151/M6/laza/ATP6AP2/beat-Va/beat-Vb/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/AP-1sigma/XNP/CG6330/beat-VI/DIP-gamma/RhoGAP100F/Vav/Ephrin/Alk/robo3/comm2/dlp/SCAR/CG30456/Unc-115a/Trim9/CG32100/gogo/shep/beat-IIIb/dpr11/Mical/pain/CG34353/lap/vvl/l(3)L1231/iav/psidin/kra/beta-Spec/meigo/mmy/kuz/Mp/Pdfr/Dscam3/conv/haf/NT1/wb/pdm3/Bsg/mid/Ptpmeg/CadN2/Src64B/Acsl/hts/Hmgcr/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/beat-IIIa/rdgC/lov/Sec15/Ten-a/SelR/ap/br/InR/Sac1/Sema5c/trol 200
GO:0006935 chemotaxis 169/3620 279/8880 0 0 0 brm/bsk/chic/eya/Dab/Lar/E(Pc)/Fas2/Fas3/FER/Glt/grn/if/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/sm/tup/Pp1-87B/ct/Ptp99A/Ptp10D/Ten-m/mew/pros/zfh1/hh/fru/otk/tok/eIF2beta/robo1/nvy/dac/Con/comm/Rac1/Cdc42/Sdc/trn/Wnt4/tutl/kermit/Sema1a/Moe/Tig/beat-Ia/shot/Cdk5/Rac2/Rho1/cher/gcm/mirr/drl/tap/CadN/Lis-1/NetA/NetB/Rab6/Vang/side/wun/fz2/14-3-3epsilon/ko/dom/Nrk/Dad/Pka-R2/caps/trio/stan/DAAM/Bap55/Cdk5alpha/jbug/bdl/Rich/beat-Ic/beat-Ib/nerfin-1/SoxN/raptor/pod1/hiw/CG6867/wgn/daw/Oli/Dscam1/Mmp2/Drl-2/ckn/Hr51/unc-5/Toll-7/Lrt/spz5/RhoGEF64C/beat-Va/beat-Vb/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/AP-1sigma/XNP/beat-VI/DIP-gamma/RhoGAP100F/Vav/Ephrin/Alk/robo3/comm2/dlp/SCAR/CG30456/Unc-115a/Trim9/gogo/beat-IIIb/dpr11/Mical/vvl/psidin/kra/beta-Spec/meigo/mmy/kuz/Mp/Dscam3/conv/haf/NT1/wb/pdm3/Bsg/mid/Ptpmeg/CadN2/Src64B/Acsl/hts/Hmgcr/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/beat-IIIa/Sec15/Ten-a/SelR/ap/InR/Sac1/Sema5c/trol 169
GO:0097485 neuron projection guidance 164/3620 270/8880 0 0 0 brm/bsk/chic/eya/Dab/Lar/E(Pc)/Fas2/Fas3/FER/Glt/grn/if/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/sm/tup/Pp1-87B/ct/Ptp99A/Ptp10D/Ten-m/mew/pros/zfh1/fru/otk/tok/eIF2beta/robo1/nvy/dac/comm/Rac1/Cdc42/Sdc/trn/Wnt4/tutl/kermit/Sema1a/Moe/Tig/beat-Ia/shot/Cdk5/Rac2/Rho1/cher/gcm/mirr/drl/tap/CadN/Lis-1/NetA/NetB/Rab6/Vang/side/fz2/14-3-3epsilon/ko/dom/Nrk/Dad/Pka-R2/caps/trio/stan/DAAM/Bap55/Cdk5alpha/jbug/bdl/Rich/beat-Ic/beat-Ib/nerfin-1/SoxN/raptor/pod1/hiw/CG6867/wgn/daw/Oli/Dscam1/Mmp2/Drl-2/ckn/Hr51/unc-5/Toll-7/spz5/RhoGEF64C/beat-Va/beat-Vb/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/AP-1sigma/XNP/beat-VI/DIP-gamma/RhoGAP100F/Vav/Ephrin/Alk/robo3/comm2/dlp/SCAR/CG30456/Unc-115a/Trim9/gogo/beat-IIIb/dpr11/Mical/vvl/psidin/kra/beta-Spec/meigo/mmy/kuz/Mp/Dscam3/conv/haf/NT1/wb/pdm3/Bsg/mid/Ptpmeg/CadN2/Src64B/Acsl/hts/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/beat-IIIa/Sec15/Ten-a/SelR/ap/InR/Sac1/Sema5c/trol 164
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "11 enrichments found"
Code
nrow(up_ekegg)
[1] 11
Code
write.csv(up_ekegg, file = paste0(output_dir, "W4KO_vs_ctrl_up_kegg.csv"))
Code
head(up_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme03013 dme03013 Nucleocytoplasmic transport - Drosophila melanogaster (fruit fly) 59/675 99/3488 0.00e+00 0.0000000 0.0000000 mago/RanGAP/Nup214/Ref1/Mtor/Hel25E/Ran/emb/Nup154/Cse1/Cbp80/Cbp20/Bin1/Dbp80/Kap-alpha1/Tnpo/mbo/msk/Kap-alpha3/Nup93-1/Nxt1/Aladin/Nup35/Ranbp21/Elys/tho2/Tnpo-SR/Nup50/Gle1/Nup54/Nup62/Rae1/Pym/thoc5/thoc7/thoc6/Srrm1/Hpr1/RnpS1/Pnn/Ranbp9/Nup93-2/Nup133/Nup98-96/Nup37/Nup358/CG18787/Arts/Apl/Ranbp16/Nup153/Snup/cdm/Nup160/Nxf3/CG34334/smt3/ebo/Gp210 59
dme03082 dme03082 ATP-dependent chromatin remodeling - Drosophila melanogaster (fruit fly) 37/675 81/3488 0.00e+00 0.0000026 0.0000023 His2Av/pho/mor/Uch-L5/Iswi/Snr1/HDAC1/Bap60/BCL7-like/CG8677/MRG15/MBD-like/Bap111/CDK2AP1/Arp8/CG31917/CG7154/YL-1/d4/CG10395/MrgBP/CG18004/DMAP1/Eaf6/polybromo/rept/pont/CG12659/Bap170/Chrac-16/Chrac-14/CG33695/Ino80/e(y)3/Mi-2/Rcd5/Caf1-55 37
dme03022 dme03022 Basal transcription factors - Drosophila melanogaster (fruit fly) 21/675 40/3488 2.50e-06 0.0000943 0.0000827 e(y)1/hay/TfIIB/TfIIFalpha/Trf/Taf1/TfIIFbeta/TfIIA-L/Taf12/Taf11/Taf2/TfIIA-S/Taf8/Taf7/Taf10b/Tfb4/Ssl1/Trf2/Xpd/Cdk7/Tfb5 21
dme03040 dme03040 Spliceosome - Drosophila melanogaster (fruit fly) 48/675 141/3488 1.88e-05 0.0005355 0.0004698 crn/SmF/l(1)10Bb/mago/Bx42/Ref1/Hel25E/noi/Spx/U2af38/ytr/Cbp80/Cbp20/CG2685/l(1)G0007/Prp18/CG4119/x16/hfp/Tango4/CG1622/CG6227/Usp39/tho2/CG3542/Phf5a/CG12343/Cypl/CG2021/CG11586/CG6610/Sf3b6/CG6418/U4-U6-60K/CG14641/Hpr1/BCAS2/c12.1/LSm3/pea/Prp19/SmD2/SmE/SNRPG/SmD1/SmB/Cdc5/caz 48
dme03250 dme03250 Viral life cycle - HIV-1 - Drosophila melanogaster (fruit fly) 17/675 33/3488 3.20e-05 0.0007289 0.0006394 BicD/Fur2/Snr1/Cdk9/Ran/Gcn5/emb/ear/fal/Hrs/Tnpo-SR/Cpsf6/Nup358/Jasper/Nup153/ALiX/nej 17
dme04330 dme04330 Notch signaling pathway - Drosophila melanogaster (fruit fly) 14/675 25/3488 4.85e-05 0.0009213 0.0008082 dsh/gro/H/Su(H)/Bx42/HDAC1/Gcn5/Atx-1/aph-1/CG6843/Nct/pen-2/nej/Psn 14
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "3 enrichments found"
Code
nrow(down_ekegg)
[1] 3
Code
write.csv(down_ekegg, file = paste0(output_dir, "W4KO_vs_ctrl_down_kegg.csv"))
Code
head(down_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme04391 dme04391 Hippo signaling pathway - fly - Drosophila melanogaster (fruit fly) 35/1081 56/3488 0.0000010 0.0001131 0.0001000 Act5C/baz/bsk/fj/ft/hth/dlg1/l(2)gl/sd/tsh/vn/mts/ex/14-3-3zeta/CycE/Gug/wts/14-3-3epsilon/jub/wgn/lft/Dlg5/mats/p53/fred/Krn/Mer/app/sdt/d/Myc/scrib/dally/wg/ds 35
dme04013 dme04013 MAPK signaling pathway - fly - Drosophila melanogaster (fruit fly) 48/1081 93/3488 0.0000216 0.0012871 0.0011385 bi/boss/bsk/chic/Jra/kay/Sos/peb/Raf/pnt/raw/sev/Btk29A/Egfr/ttk/vap/vn/mts/aos/pros/drk/14-3-3zeta/Dsor1/Rac1/msn/eff/Mef2/Ask1/Rac2/RhoL/sty/kek1/ksr/p38a/PTP-ER/edl/Mkk4/Mekk1/Tak1/slpr/wgn/CYLD/Krn/Tab2/alph/puc/cic/Src64B 48
dme04080 dme04080 Neuroactive ligand-receptor interaction - Drosophila melanogaster (fruit fly) 33/1081 67/3488 0.0012238 0.0485425 0.0429390 mAChR-A/baz/alphaTry/5-HT1A/Rdl/5-HT7/GluRIA/TkR99D/TkR86C/Lcch3/deltaTry/Nmdar1/zetaTry/Dop1R1/mGluR/GABA-B-R2/MSBP/CCKLR-17D3/GABA-B-R3/CCHa2-R/Dh44-R2/Dh44-R1/CrzR/mAChR-B/CCAP-R/lambdaTry/kappaTry/Dh31-R/Nmdar2/Dop2R/GABA-B-R1/5-HT1B/GluRIB 33

3.6 O1_KO

Code
results_O1KO_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "O1_KO", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_O1_KO_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_O1KO_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 2827
Code
### MA plot
generate_ma(res)
Warning: Removed 1 rows containing missing values (`geom_point()`).

Code
### Volcano plot
generate_volcano(res, "O1_KO", "Ctrl")
Warning: Removed 1 rows containing missing values (`geom_point()`).

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 1027
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "199 enrichments found"
Code
nrow(ego)
[1] 199
Code
write.csv(ego, file = paste0(output_dir, "O1KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0006259 DNA metabolic process 96/841 329/8880 0 0 0 aurA/fs(1)Ya/gnu/hay/Hmr/lds/mei-9/mei-218/mago/mu2/mus101/DNApol-zeta/mus304/mus312/otu/phr/plu/rec/spn-B/spn-D/spn-E/Su(var)205/Cdk1/tam/piwi/PCNA/CycE/Mlh1/Pms2/RnrL/Irbp/mof/Orc5/spel1/tos/lok/mre11/rnh1/Orc6/CG14803/CG14814/Drep1/Rad9/CG6171/Cdc7/Ercc1/mei-217/CG4078/Caf1-180/CG2990/Psf3/Arp8/CG13690/HP6/Nse1/Mnn1/RfC3/RPA2/CG8679/mms4/CG9062/Ctf4/nopo/CG11788/Mes4/rad50/Psf1/Sse/CG3448/Msh6/HipHop/Snm1/dgrn/CG10445/CG9601/Fancl/BigH1/Fancd2/CG6985/Top3alpha/RecQ4/mus81/Ku80/Chrac-14/tefu/Tengl3/CG31807/SMC5/hd/gkt/RfC4/Xpd/Psf2/ver/Caf1-55/RPA3 96
GO:0006974 cellular response to DNA damage stimulus 74/841 233/8880 0 0 0 aurA/hay/lds/mei-9/mei-218/mu2/mus101/DNApol-zeta/mus304/mus312/phr/spn-B/spn-D/Cdk1/tam/PCNA/Dsor1/Mlh1/Pms2/Irbp/l(2)dtl/mof/spel1/tos/lok/mre11/CG14803/CG14814/Rad9/CG6171/Cdc7/Ercc1/CG4078/Corp/Caf1-180/CG8173/Arp8/CG13690/Nse1/Mnn1/RfC3/RPA2/CG8679/mms4/CG9062/Ctf4/rad50/CG3448/Msh6/HipHop/Snm1/CG1218/dgrn/CG10445/CG9601/Fancl/Fancd2/CG6985/Top3alpha/RecQ4/mus81/armi/Ku80/Chrac-14/tefu/CG31807/Claspin/SMC5/gkt/RfC4/grp/Xpd/Psf2/RPA3 74
GO:0006281 DNA repair 62/841 186/8880 0 0 0 aurA/hay/lds/mei-9/mei-218/mu2/DNApol-zeta/mus304/mus312/phr/spn-B/spn-D/tam/PCNA/Mlh1/Pms2/Irbp/mof/spel1/tos/lok/mre11/CG14803/CG14814/Rad9/CG6171/Cdc7/Ercc1/CG4078/Caf1-180/Arp8/CG13690/Nse1/Mnn1/RfC3/RPA2/CG8679/mms4/CG9062/Ctf4/rad50/CG3448/Msh6/Snm1/dgrn/CG10445/CG9601/Fancl/Fancd2/CG6985/Top3alpha/RecQ4/mus81/Ku80/tefu/CG31807/SMC5/gkt/RfC4/Xpd/Psf2/RPA3 62
GO:0051276 chromosome organization 79/841 298/8880 0 0 0 aub/CycA/CycB/fs(1)Ya/Hmr/mei-218/mu2/mus301/mus304/rec/spn-B/spn-D/spn-E/Su(var)3-7/Su(var)205/Klp67A/shtd/rhi/wapl/piwi/gammaTub37C/Irbp/prod/Orc5/spel1/pgc/lok/mre11/san/CG11403/ball/CG4078/APC7/Fign/HP6/lid/RfC3/Cdc23/Nup50/Mos/tej/krimp/CG11788/Panx/Rae1/rad50/Sse/mad2/PCID2/HipHop/Snm1/Hat1/CG8223/BigH1/cal1/Fancd2/Nup133/Mink/Top3alpha/RecQ4/ida/Ku80/Chrac-14/SA-2/tefu/Haspin/Nnf1b/Zwilch/del/gwl/gammaTub23C/cuff/RfC4/Incenp/Xpd/Psf2/ver/BubR1/vas 79
GO:1903046 meiotic cell cycle process 61/841 201/8880 0 0 0 aub/Axs/cup/Grip91/lds/twe/mei-9/mei-218/mus301/mus304/mus312/rec/RanGAP/spn-B/spn-D/spn-E/Cdk1/Klp67A/yem/gammaTub37C/Cks30A/mtrm/dhd/Pp2B-14D/mael/lok/Ran/mre11/Grip75/Grip128/ball/Ercc1/mei-217/HP6/lid/Grip71/CG8679/mms4/Mos/tej/krimp/Rae1/rad50/Sse/Syx13/Msh6/cal1/Fancd2/CG6985/mus81/tefu/del/tef/sra/gwl/gammaTub23C/wisp/cuff/Incenp/BubR1/vas 61
GO:0051321 meiotic cell cycle 62/841 215/8880 0 0 0 aub/Axs/cup/Grip91/lds/twe/mei-9/mei-218/mus301/mus304/mus312/rec/RanGAP/spn-B/spn-D/spn-E/Cdk1/Klp67A/yem/gammaTub37C/Cks30A/mtrm/dhd/Pp2B-14D/mael/lok/Ran/mre11/mei-P26/Grip75/Grip128/ball/Ercc1/mei-217/HP6/lid/Grip71/CG8679/mms4/Mos/tej/krimp/Rae1/rad50/Sse/Syx13/Msh6/cal1/Fancd2/CG6985/mus81/tefu/del/tef/sra/gwl/gammaTub23C/wisp/cuff/Incenp/BubR1/vas 62
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 4169
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "75 enrichments found"
Code
nrow(ego)
[1] 75
Code
write.csv(ego, file = paste0(output_dir, "O1KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0040003 chitin-based cuticle development 107/3366 147/8880 0e+00 0.0000000 0.0000000 amd/ck/dib/ecd/Gld/kkv/knk/l(2)gl/Lcp4/Pu/Ccp84Ag/Ccp84Ae/Ccp84Ad/Ccp84Ac/Ccp84Ab/Ccp84Aa/Acp1/fz2/Idgf1/pio/Lcp65Ag2/Lcp65Ag2/Lcp65Af/Lcp65Ac/Lcp65Ab2/Alas/Cht2/Gasp/Axn/Vajk3/Cpr35B/Vajk1/gammaCOP/Chs2/Cpr11A/obst-A/obst-E/TwdlE/Vajk2/Hakai/Cpr47Ef/Twdlbeta/Cpr49Ac/Cpr49Ae/Cpr49Ag/Cpr50Cb/Cpr51A/Cpr56F/Cpr62Bb/Cpr62Bc/Cht7/Cpr64Ab/Cpr64Ad/Cpr65Ea/Cpr65Ec/Cpr66Cb/Cpr67B/Cpr67Fa1/Cpr67Fa2/Cpr76Bc/Cpr78Cb/Cpr78Cc/TwdlG/Cht5/TwdlW/Cpr92F/TwdlM/TwdlP/TwdlB/TwdlL/TwdlN/TwdlD/TwdlQ/TwdlC/Cpr97Ea/Cpr97Eb/Cpr100A/Cpr65Ax2/Vajk4/TwdlH/fred/Cpr66D/Cpr65Av/TwdlY/TwdlX/Twdlalpha/dpy/Cpr31A/Idgf5/pan/Cpr47Eg/Lcp65Ag3/puc/Tb/CG34461/pot/mmy/drd/serp/wol/mgl/verm/dsx-c73A/sdt/Sec24CD/Eogt/stw 107
GO:0042335 cuticle development 126/3366 185/8880 0e+00 0.0000000 0.0000000 amd/ck/Ddc/dib/ect/ecd/Gld/kkv/knk/l(2)gl/Lcp4/m/Pu/sha/Ser/Ccp84Ag/Ccp84Ae/Ccp84Ad/Ccp84Ac/Ccp84Ab/Ccp84Aa/ple/Cyp4g1/fng/Acp1/p38a/fz2/Idgf1/pio/Lcp65Ag2/Lcp65Ag2/Lcp65Af/Lcp65Ac/Lcp65Ab2/Alas/Cht2/WASp/Mkk4/Gasp/gus/Axn/Vajk3/Cpr35B/Vajk1/gammaCOP/Chs2/Cpr11A/obst-A/snsl/obst-E/TwdlE/Vajk2/Hakai/Cpr47Ef/Twdlbeta/Cpr49Ac/Cpr49Ae/Cpr49Ag/Cpr50Cb/Cpr51A/Cpr56F/Cpr62Bb/Cpr62Bc/Cht7/Cpr64Ab/Cpr64Ad/Cpr65Ea/Cpr65Ec/Cpr66Cb/Cpr67B/Cpr67Fa1/Cpr67Fa2/wls/Cpr76Bc/Cpr78Cb/Cpr78Cc/TwdlG/Cht5/TwdlW/Cpr92F/TwdlM/TwdlP/TwdlB/TwdlL/TwdlN/TwdlD/TwdlQ/TwdlC/Cpr97Ea/Cpr97Eb/snu/Cpr100A/boi/Cpr65Ax2/FASN2/Vajk4/TwdlH/fred/Cpr66D/Cpr65Av/TwdlY/TwdlX/Twdlalpha/dpy/Cpr31A/Idgf5/pan/Cpr47Eg/Lcp65Ag3/Desat1/puc/Tb/CG34461/pot/mmy/drd/serp/wol/mgl/verm/dsx-c73A/sdt/f/Sec24CD/Eogt/stw 126
GO:0008045 motor neuron axon guidance 51/3366 75/8880 1e-07 0.0001060 0.0001033 dsh/Fas2/for/Glt/grn/Lim3/spz/tup/Ten-m/zfh1/robo1/Sdc/Wnt4/fra/beat-Ia/Cdk5/Ptp69D/Rac2/Rho1/cher/NetA/NetB/side/fz2/ko/caps/Cdk5alpha/beat-Ic/fend/CG6867/beat-IIIc/Oli/Mmp2/Drl-2/unc-5/spz5/beat-Va/beat-IIb/beat-IIa/beat-VI/beat-IIIb/vvl/Mp/CG42637/haf/NT1/dally/Nrg/beat-IIIa/CG42637/trol 51
GO:0099536 synaptic signaling 155/3366 296/8880 2e-07 0.0001117 0.0001088 Ace/nAChRalpha1/mAChR-A/nAChRbeta1/nAChRalpha2/boss/Cat/ChAT/comt/dnc/Fas2/Galphas/mr/na/Atpalpha/pum/rt/rut/sesB/Sh/slo/svp/nAChRbeta2/5-HT1A/Csp/Syt1/Rdl/Ten-m/5-HT7/Rop/Syn/Cbp53E/TkR99D/Rap1/Rab3/CanB/Lcch3/Nmdar1/Dcp-1/Snap25/Dop1R1/nSyb/stnB/mGluR/dtr/wit/b6/Oamb/GluClalpha/Srpk79D/X11L/Rep/Eaat1/Rab3-GAP/GABA-B-R2/Syt4/GluRIID/nAChRalpha5/NtR/Vap33/Pink1/BTBD9/hiw/Rab3-GEF/Frq1/wgn/VGlut/Nlg2/CG7781/CG14275/nAChRalpha6/atilla/CG9336/CG9338/Dyb/CG6329/VGAT/Dg/DAT/P32/CG10339/spz5/unc-13-4A/PGRP-LC/Smn/rogdi/Neurochondrin/mAChR-B/Kcmf1/DCAF12/KaiR1D/Nrx-1/spas/Tusp/CDase/kek6/Syt7/Ekar/CG11155/Scamp/Alk/cpx/Fsn/GluRIIC/RIC-3/Nlg1/GluRIIE/CG31475/CG31675/CG31676/scramb1/Csas/CG32447/be/X11Lbeta/Rim/Cadps/Frq2/Nlg3/Nlg4/stmA/beta-Spec/alpha-Spec/Pka-R1/brp/Ca-beta/Hip14/Dys/qvr/Atg1/stj/Syt12/Sytalpha/Sytbeta/NT1/Exn/nemy/Ank2/Rbp/Shab/cac/5-HT1B/Acsl/Ih/nwk/cpo/dysc/Slob/Fife/Neto/nAChRalpha4/Sec8/Esyt2/VAChT/cv-c 155
GO:0099537 trans-synaptic signaling 153/3366 292/8880 2e-07 0.0001117 0.0001088 Ace/nAChRalpha1/mAChR-A/nAChRbeta1/nAChRalpha2/boss/Cat/ChAT/comt/dnc/Fas2/Galphas/mr/na/Atpalpha/pum/rt/rut/sesB/Sh/slo/svp/nAChRbeta2/5-HT1A/Csp/Syt1/Rdl/Ten-m/5-HT7/Rop/Syn/Cbp53E/TkR99D/Rap1/Rab3/CanB/Lcch3/Nmdar1/Dcp-1/Snap25/Dop1R1/nSyb/stnB/mGluR/dtr/wit/b6/Oamb/GluClalpha/Srpk79D/X11L/Rep/Eaat1/Rab3-GAP/GABA-B-R2/Syt4/GluRIID/nAChRalpha5/NtR/Vap33/Pink1/BTBD9/hiw/Rab3-GEF/Frq1/wgn/VGlut/Nlg2/CG7781/CG14275/nAChRalpha6/atilla/CG9336/CG9338/Dyb/CG6329/VGAT/Dg/DAT/P32/CG10339/spz5/unc-13-4A/PGRP-LC/Smn/rogdi/Neurochondrin/mAChR-B/DCAF12/KaiR1D/Nrx-1/spas/CDase/kek6/Syt7/Ekar/CG11155/Scamp/Alk/cpx/Fsn/GluRIIC/RIC-3/Nlg1/GluRIIE/CG31475/CG31675/CG31676/scramb1/Csas/CG32447/be/X11Lbeta/Rim/Cadps/Frq2/Nlg3/Nlg4/stmA/beta-Spec/alpha-Spec/Pka-R1/brp/Ca-beta/Hip14/Dys/qvr/Atg1/stj/Syt12/Sytalpha/Sytbeta/NT1/Exn/nemy/Ank2/Rbp/Shab/cac/5-HT1B/Acsl/Ih/nwk/cpo/dysc/Slob/Fife/Neto/nAChRalpha4/Sec8/Esyt2/VAChT/cv-c 153
GO:0003012 muscle system process 27/3366 33/8880 3e-07 0.0001117 0.0001088 Ca-alpha1D/Mlc1/Mlc2/Tm2/up/RyR/Ms/CG3078/Eaat1/Amph/CG4928/slmo/KCNQ/Chd64/Neurochondrin/CG6972/robo3/Unc-89/Lasp/Nox/Mp/Mhc/tn/Fhos/SclA/SclB/wupA 27
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "12 enrichments found"
Code
nrow(up_ekegg)
[1] 12
Code
write.csv(up_ekegg, file = paste0(output_dir, "O1KO_vs_ctrl_up_kegg.csv"))
Code
head(up_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme03460 dme03460 Fanconi anemia pathway - Drosophila melanogaster (fruit fly) 13/294 26/3488 0.0000000 0.0000034 0.0000031 mei-9/DNApol-zeta/spn-D/Mlh1/Pms2/Ercc1/Usp1/mms4/CG9062/Fancl/Fancd2/Top3alpha/mus81 13
dme03050 dme03050 Proteasome - Drosophila melanogaster (fruit fly) 18/294 53/3488 0.0000001 0.0000063 0.0000057 Prosbeta6/Prosalpha4/Prosbeta1/Prosalpha5/Prosalpha7/Prosbeta3/Rpn12/Prosbeta5/Prosbeta4/Pomp/Rpt3R/Rpt6R/CG30382/Prosalpha2/Prosbeta7/Prosalpha6/Prosalpha3/Prosalpha1 18
dme03440 dme03440 Homologous recombination - Drosophila melanogaster (fruit fly) 9/294 23/3488 0.0000539 0.0019209 0.0017385 mus101/spn-B/spn-D/mre11/mms4/rad50/Top3alpha/mus81/tefu 9
dme03430 dme03430 Mismatch repair - Drosophila melanogaster (fruit fly) 8/294 19/3488 0.0000768 0.0020547 0.0018596 PCNA/Mlh1/Pms2/spel1/tos/RfC3/Msh6/RfC4 8
dme03013 dme03013 Nucleocytoplasmic transport - Drosophila melanogaster (fruit fly) 19/294 99/3488 0.0004581 0.0098036 0.0088729 mago/RanGAP/Ref1/Ran/Cse1/Nup35/Tnpo-SR/Nup50/Nup54/Rae1/Pym/thoc6/Ranbp9/Nup133/Ndc1/Nup37/Apl/Nxf3/smt3 19
dme00240 dme00240 Pyrimidine metabolism - Drosophila melanogaster (fruit fly) 10/294 38/3488 0.0008693 0.0155022 0.0140305 awd/Dhod/RnrL/dnk/Ts/nmdyn-D6/CG5757/CG3788/cN-IIIB/dUTPase 10
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "9 enrichments found"
Code
nrow(down_ekegg)
[1] 9
Code
write.csv(down_ekegg, file = paste0(output_dir, "O1KO_vs_ctrl_down_kegg.csv"))
Code
head(down_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme04512 dme04512 ECM-receptor interaction - Drosophila melanogaster (fruit fly) 10/1011 11/3488 0.0000329 0.0039189 0.0035012 Col4a1/LanA/vkg/Hml/CG3168/Tsp/Dg/wb/LanB1/LanB2 10
dme04080 dme04080 Neuroactive ligand-receptor interaction - Drosophila melanogaster (fruit fly) 34/1011 67/3488 0.0001256 0.0074760 0.0066791 mAChR-A/baz/alphaTry/5-HT1A/Rdl/5-HT7/GluRIA/TkR99D/Lcch3/betaTry/deltaTry/gammaTry/Nmdar1/epsilonTry/etaTry/zetaTry/Dop1R1/Try29F/mGluR/GABA-B-R2/CCKLR-17D3/GABA-B-R3/Dh44-R2/CrzR/mAChR-B/lambdaTry/kappaTry/Dh31-R/Nmdar2/Dop2R/Gpb5/GABA-B-R1/5-HT1B/GluRIB 34
dme04330 dme04330 Notch signaling pathway - Drosophila melanogaster (fruit fly) 16/1011 25/3488 0.0002794 0.0110812 0.0099001 Dl/dsh/dx/gro/H/mam/Ser/Su(H)/fng/HDAC1/Gcn5/Hey/aph-1/CG6843/Tace/Psn 16
dme03013 dme03013 Nucleocytoplasmic transport - Drosophila melanogaster (fruit fly) 43/1011 99/3488 0.0013278 0.0395031 0.0352925 eEF1alpha2/sbr/Nup214/Mtor/Hel25E/emb/Nup154/Cbp80/Cbp20/Bin1/Kap-alpha1/Tnpo/msk/Nup93-1/Nxt1/Ranbp21/Elys/tho2/Nup44A/Gle1/tsu/Phax/Nup188/Nup62/Nup75/thoc7/Srrm1/Hpr1/RnpS1/Pnn/CG14712/Nup93-2/Nup43/Nup58/Nup98-96/Nup358/Ranbp16/Snup/Nup160/Fs(2)Ket/Kr-h2/ebo/Gp210 43
dme00230 dme00230 Purine metabolism - Drosophila melanogaster (fruit fly) 40/1011 92/3488 0.0018919 0.0450271 0.0402276 Gart/Aprt/dnc/Pgm1/pn/rut/ry/Ac76E/Gyc32E/RnrS/Gycalpha99B/Gycbeta100B/Papss/Paics/Adk1/Ac3/NTPase/AdSS/nmdyn-D7/CG6106/CG8891/CG11050/Pgm2a/veil/CG16758/Adgf-A/Pde6/Gyc88E/Gyc89Db/AdSL/ACXD/NT5E-2/AMPdeam/CG34357/Pde9/CG42637/DhpD/Pde1c/Gyc76C/Pde8 40
dme03022 dme03022 Basal transcription factors - Drosophila melanogaster (fruit fly) 20/1011 40/3488 0.0038833 0.0664614 0.0593773 e(y)1/Taf4/Trf/Taf1/Taf5/TfIIA-L/Taf12/Taf2/TfIIEalpha/CycH/Taf7/Mat1/bip2/Taf10/Tfb4/Tfb1/Ssl1/mrn/Trf2/Tfb5 20

3.7 O1W4_KO

Code
results_O1W4KO_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "O1W4_KO", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_O1W4_KO_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_O1W4KO_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 4881
Code
### MA plot
generate_ma(res)

Code
### Volcano plot
generate_volcano(res, "O1W4_KO", "Ctrl")

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 2514
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "313 enrichments found"
Code
nrow(ego)
[1] 313
Code
write.csv(ego, file = paste0(output_dir, "O1W4KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0006397 mRNA processing 136/2093 308/8880 0 0 0 bru1/aub/crn/SmF/ecd/l(1)10Bb/mago/mod/nos/snf/Bx42/U2af50/Rox8/Ref1/Gem3/Hel25E/noi/Clp/Spx/U2af38/ytr/ubl/Cbp80/Cbp20/sip1/Srp54/Bin1/Neos/Cpsf160/CG2685/SF1/l(1)G0007/CG5808/Prp18/CstF64/Cpsf100/Nab2/CG4119/x16/hfp/spt4/CG16865/icln/asf1/CG1677/Tango4/CG1622/CG5877/CG6227/CG9213/CG6179/Usp39/cactin/tho2/Cwc25/Tnpo-SR/CG3542/Phf5a/Mettl14/Bug22/CG6686/CG10466/Ars2/PAN2/tsu/Lsm11/Lsm10/CG12343/FLASH/cbc/CG11808/CG8435/CG7747/CG15084/CG4266/thoc7/Cypl/CG7974/CG2021/CG13807/PAN3/CG11586/CG6610/Sf3b6/ldbr/Cpsf6/CG3335/CG10907/CG10418/Sf3a2/Srrm1/Smn/U4-U6-60K/Gem2/CG14641/Fip1/CG1965/Tailor/CG7483/beag/RnpS1/Pnn/CG10324/mdlc/tbrd-1/PQBP1/CG6066/BCAS2/Moca-cyp/CstF50/CG17454/c12.1/CG5986/Prp38/LSm3/snama/Spf45/ncm/pea/rig/Larp7/wisp/Rbp1/Prp19/ymp/SmD2/SmE/SNRPG/snRNP-U1-C/SmD1/SmB/Rpb4/Cdc5/SF2/caz/Rnmt 136
GO:0051276 chromosome organization 132/2093 298/8880 0 0 0 aub/c(3)G/Cam/CycA/CycB/Hmr/kel/pho/mei-S332/mod(mdg4)/mr/ms(3)K81/mu2/mus301/mus304/Blm/ncd/pim/polo/spn-B/spn-D/spn-E/sub/Su(var)3-7/Su(var)205/Su(var)2-10/thr/Klp61F/Klp67A/shtd/rhi/wapl/piwi/Top1/gammaTub37C/Parp/eff/Mtor/prod/glu/Cul1/Lis-1/HDAC1/pgc/Stat92E/mre11/Mcm7/Nup154/DCTN2-p50/san/aurB/SuUR/Bub3/Bap60/SkpA/Cdc16/Mcm6/cav/Mad1/Nipped-B/Hus1-like/Eb1/Kmn1/MRG15/Nup93-1/SMC2/ball/MBD-like/RfC38/APC7/Fign/Pif1/HP6/lid/Nuf2/baf/borr/RfC3/Mcm10/CG15237/Nup50/Mos/mars/tej/krimp/Nup62/Oxp/Panx/Rae1/mad2/nudE/PCID2/Smn/Snm1/Hat1/Arl8/CG8223/BigH1/Mau2/cal1/Nup133/HP1c/Nup98-96/Mink/Kat60/cana/Chrac-14/SA-2/Haspin/Rif1/nenya/Nnf1b/His1:CG33801/His1:CG33861/tum/egg/Spc25/lmgA/Su(var)3-3/gwl/vtd/Incenp/zuc/CoRest/Psf2/Mi-2/ver/BubR1/Topors/Yeti/solo/Top2 132
GO:0016071 mRNA metabolic process 160/2093 395/8880 0 0 0 bru1/aub/crn/SmF/ecd/l(1)10Bb/mago/mod/nos/pum/snf/spn-E/me31B/RpII15/Bx42/piwi/TfIIB/U2af50/Rox8/Ref1/Gem3/Hel25E/noi/Clp/Spx/smg/U2af38/ytr/ubl/Cbp80/Cbp20/sip1/Srp54/Bin1/Neos/Cpsf160/CG2685/SF1/l(1)G0007/CG5808/Prp18/CstF64/Cpsf100/Nab2/CG4119/x16/hfp/spt4/CG16865/icln/asf1/CG1677/Tango4/CG1622/CG5877/CG6227/CG9213/CG6179/Usp39/cactin/tho2/Cwc25/Tnpo-SR/CG3542/Phf5a/Rat1/Mettl14/Bug22/CG6686/CG10466/Ars2/PAN2/tsu/Lsm11/Lsm10/CG12343/FLASH/cbc/CG11808/CG8435/CG7747/CG15084/CG4266/CG3045/Pym/DCP1/Upf3/thoc7/Cypl/CG7974/CG2021/CG13807/PAN3/CG11586/CG6610/Sf3b6/ldbr/Cpsf6/CG3335/CG10907/CG10418/Sf3a2/Srrm1/Smn/CG9705/U4-U6-60K/Edc3/Gem2/CG14641/Fip1/CG2091/CG1965/Tailor/CG7483/beag/RnpS1/Pnn/CG10324/mdlc/tbrd-1/Ythdf/PQBP1/CG6066/BCAS2/Moca-cyp/CstF50/Marf1/CG17454/c12.1/armi/CG5986/Prp38/Rpb7/LSm3/CG31643/CG31759/gw/snama/Spf45/ncm/pea/AGO3/rig/Larp7/wisp/Rbp1/Prp19/ymp/SmD2/SmE/SNRPG/snRNP-U1-C/SmD1/SmB/Rpb4/Cdc5/Patr-1/SF2/caz/Rnmt 160
GO:0008380 RNA splicing 120/2093 270/8880 0 0 0 bru1/aub/crn/SmF/ecd/l(1)10Bb/mago/mod/snf/Bx42/U2af50/Rox8/Parp/Ref1/Gem3/Hel25E/noi/Spx/U2af38/ytr/ubl/Cbp80/Cbp20/Parg/sip1/Srp54/Bin1/Neos/ssx/CG2685/SF1/l(1)G0007/CG5808/Prp18/CG4119/x16/hfp/CG16865/icln/asf1/CG1677/Tango4/CG1622/CG5877/CG6227/CG9213/CG6179/Usp39/cactin/Cwc25/Tnpo-SR/CG3542/Phf5a/Mettl14/Bug22/CG6686/CG10466/Ars2/tsu/Lsm11/CG12343/cbc/CG11808/CG8435/CG7747/CG15084/Cypl/CG7974/CG2021/CG11586/CG6610/Sf3b6/ldbr/BI-1/Cpsf6/CG3335/CG10907/Tsen54/CG10418/Sf3a2/Srrm1/Smn/U4-U6-60K/Gem2/CG14641/CG1965/CG7483/beag/RnpS1/Pnn/CG10324/mdlc/tbrd-1/PQBP1/CG6066/BCAS2/Moca-cyp/CG17454/c12.1/CG5986/Prp38/LSm3/Tsen34/Spf45/ncm/pea/rig/Larp7/Rbp1/Prp19/ymp/SmD2/SmE/SNRPG/snRNP-U1-C/SmD1/SmB/Cdc5/SF2/caz 120
GO:0000375 RNA splicing, via transesterification reactions 113/2093 250/8880 0 0 0 bru1/aub/crn/SmF/ecd/l(1)10Bb/mago/mod/snf/Bx42/U2af50/Rox8/Ref1/Gem3/Hel25E/noi/Spx/U2af38/ytr/ubl/Cbp80/Cbp20/sip1/Srp54/Bin1/Neos/CG2685/SF1/l(1)G0007/CG5808/Prp18/CG4119/x16/hfp/CG16865/icln/asf1/CG1677/Tango4/CG1622/CG5877/CG6227/CG9213/CG6179/Usp39/cactin/Cwc25/Tnpo-SR/CG3542/Phf5a/Mettl14/Bug22/CG6686/CG10466/Ars2/tsu/Lsm11/CG12343/CG11808/CG8435/CG7747/CG15084/Cypl/CG7974/CG2021/CG11586/CG6610/Sf3b6/ldbr/Cpsf6/CG3335/CG10907/CG10418/Sf3a2/Srrm1/Smn/U4-U6-60K/Gem2/CG14641/CG1965/CG7483/beag/RnpS1/Pnn/CG10324/mdlc/tbrd-1/PQBP1/CG6066/BCAS2/Moca-cyp/CG17454/c12.1/CG5986/Prp38/LSm3/Spf45/ncm/pea/rig/Larp7/Rbp1/Prp19/ymp/SmD2/SmE/SNRPG/snRNP-U1-C/SmD1/SmB/Cdc5/SF2/caz 113
GO:0000377 RNA splicing, via transesterification reactions with bulged adenosine as nucleophile 113/2093 250/8880 0 0 0 bru1/aub/crn/SmF/ecd/l(1)10Bb/mago/mod/snf/Bx42/U2af50/Rox8/Ref1/Gem3/Hel25E/noi/Spx/U2af38/ytr/ubl/Cbp80/Cbp20/sip1/Srp54/Bin1/Neos/CG2685/SF1/l(1)G0007/CG5808/Prp18/CG4119/x16/hfp/CG16865/icln/asf1/CG1677/Tango4/CG1622/CG5877/CG6227/CG9213/CG6179/Usp39/cactin/Cwc25/Tnpo-SR/CG3542/Phf5a/Mettl14/Bug22/CG6686/CG10466/Ars2/tsu/Lsm11/CG12343/CG11808/CG8435/CG7747/CG15084/Cypl/CG7974/CG2021/CG11586/CG6610/Sf3b6/ldbr/Cpsf6/CG3335/CG10907/CG10418/Sf3a2/Srrm1/Smn/U4-U6-60K/Gem2/CG14641/CG1965/CG7483/beag/RnpS1/Pnn/CG10324/mdlc/tbrd-1/PQBP1/CG6066/BCAS2/Moca-cyp/CG17454/c12.1/CG5986/Prp38/LSm3/Spf45/ncm/pea/rig/Larp7/Rbp1/Prp19/ymp/SmD2/SmE/SNRPG/snRNP-U1-C/SmD1/SmB/Cdc5/SF2/caz 113
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 4418
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "511 enrichments found"
Code
nrow(ego)
[1] 511
Code
write.csv(ego, file = paste0(output_dir, "O1W4KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0031175 neuron projection development 296/3669 472/8880 0 0 0 acj6/Appl/brm/bsk/chif/chic/eya/D/dnc/EcR/ems/en/E(spl)m8-HLH/Fas2/Fas3/FER/Glt/grn/if/kay/kn/Arpc1/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/Poxm/ras/rut/Scm/Sh/shn/sm/ss/Tm1/trc/ttk/ptc/tup/hid/oc/Nrt/ct/Wnt2/Ptp4E/Ptp99A/Ptp10D/Galphaq/Ten-m/mew/Rop/pros/zfh1/rho/N/fru/sim/otk/ci/tok/eIF2beta/Mbs/Sox14/robo1/nvy/dac/Con/comm/Wnt5/hep/Cdc42/Akt1/SrpRalpha/Ror/Sdc/trn/Wnt4/tutl/kermit/sqz/TrpRS/msn/Sema1a/babo/fra/mas/Med/Moe/Tig/wts/Hem/Ret/beat-Ia/shot/Cdk5/Rac2/Rel/Rho1/cher/gcm/mirr/Pi3K92E/Shc/chn/drl/tap/CadN/Dref/Nak/NetA/NetB/Rab6/Vang/ara/nompA/side/fz2/gcm2/ko/Rfx/Nrk/Dad/chb/nudC/Tor/Pka-R2/caps/RhoGEF2/wit/WASp/stan/Optix/DAAM/TBPH/wrapper/Rok/RhoGAPp190/Pten/jim/Cdk5alpha/jbug/bdl/Rich/beat-Ic/Fmr1/seq/nerfin-1/SoxN/ssh/lva/pod1/Dlic/CG11584/hiw/HUWE1/Efhc1.1/CG6867/wgn/rictor/erm/CG14043/BuGZ/Oli/Hasp/CdGAPr/Dscam1/Mmp2/Drl-2/CG6701/ckn/unc-5/Toll-7/Efhc1.2/Liprin-gamma/babos/CG13531/Mmp1/p130CAS/spz5/RhoGEF64C/frac/path/Toll-6/sff/wnd/CG8765/Cont/kat-60L1/Neurochondrin/beat-Vc/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/Dcr-1/atl/XNP/beat-VI/DIP-gamma/RhoGAP100F/mav/Vav/SMC1/Alk/robo3/exex/comm2/dlp/SCAR/bchs/ham/Liprin-alpha/CG30456/Unc-115a/Trim9/gogo/LRP1/CG33111/sav/dpr11/Mical/RunxA/nub/spri/sif/Snoo/jing/jeb/vvl/psidin/kra/beta-Spec/alpha-Spec/futsch/Mob4/crb/mmy/cwo/kuz/Mp/mfas/Dscam3/haf/NT1/wb/pdm3/Dbx/Ank2/Dhc64C/Bsg/mid/Ire1/Ptpmeg/CadN2/Sec23/dati/Src64B/cnc/Acsl/Patronin/Dscam4/dar1/hts/Dark/sick/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/Dscam2/Mtp/Kr-h1/tau/Ten-a/unc-104/SelR/Frl/ap/br/S6k/InR/Sac1/Sema5c/Adf1/trol/btl/sbb/Imp/cv-c 296
GO:0048667 cell morphogenesis involved in neuron differentiation 273/3669 430/8880 0 0 0 acj6/Appl/baz/brm/bsk/chif/chic/eya/D/dnc/EcR/ems/en/E(spl)m8-HLH/Fas2/Fas3/FER/Glt/grn/if/kay/kn/Arpc1/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/Poxm/ras/rut/Scm/Sh/shn/sm/ss/Tm1/trc/ptc/tup/hid/oc/Nrt/ct/Wnt2/Ptp4E/Ptp99A/Ptp10D/Galphaq/Ten-m/mew/pros/zfh1/rho/N/fru/sim/otk/ci/tok/eIF2beta/Sox14/robo1/nvy/dac/comm/Wnt5/hep/Cdc42/SrpRalpha/Sdc/trn/Wnt4/tutl/kermit/sqz/TrpRS/msn/Sema1a/babo/fra/Med/Moe/Tig/wts/Hem/beat-Ia/shot/Cdk5/Rac2/Rel/Rho1/cher/gcm/mirr/Pi3K92E/Shc/chn/drl/tap/CadN/Dref/Nak/NetA/NetB/Rab6/Vang/ara/nompA/side/fz2/gcm2/ko/Rfx/Nrk/Dad/chb/nudC/Tor/Pka-R2/caps/RhoGEF2/wit/WASp/stan/Optix/DAAM/Rok/RhoGAPp190/jim/Cdk5alpha/jbug/bdl/Rich/beat-Ic/Fmr1/seq/nerfin-1/SoxN/ssh/lva/pod1/Dlic/CG11584/hiw/HUWE1/Efhc1.1/CG6867/wgn/rictor/erm/CG14043/BuGZ/Oli/CdGAPr/Dscam1/Mmp2/Drl-2/CG6701/ckn/unc-5/Toll-7/Liprin-gamma/babos/CG13531/spz5/RhoGEF64C/frac/path/Toll-6/sff/wnd/CG8765/Cont/kat-60L1/beat-Vc/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/Dcr-1/atl/XNP/beat-VI/DIP-gamma/RhoGAP100F/Vav/SMC1/Alk/robo3/exex/comm2/dlp/SCAR/ham/Liprin-alpha/CG30456/Unc-115a/Trim9/gogo/LRP1/sav/dpr11/Mical/Patj/RunxA/nub/spri/sif/Snoo/jing/jeb/vvl/psidin/kra/beta-Spec/alpha-Spec/futsch/crb/mmy/cwo/kuz/Mp/mfas/Dscam3/haf/NT1/wb/pdm3/Dbx/Ank2/Dhc64C/Bsg/mid/Ire1/Ptpmeg/CadN2/Sec23/Src64B/cnc/Acsl/Dscam4/dar1/hts/Dark/sick/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/Dscam2/Mtp/tau/Ten-a/unc-104/SelR/Frl/ap/br/S6k/InR/Sac1/Sema5c/trol/btl/sbb/cv-c 273
GO:0048858 cell projection morphogenesis 276/3669 438/8880 0 0 0 acj6/Appl/brm/bsk/chif/chic/eya/D/dnc/EcR/ems/en/E(spl)m8-HLH/Fas2/Fas3/FER/Glt/grn/if/kay/kn/Arpc1/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/Poxm/ras/rut/Scm/Sh/shn/sm/ss/Tm1/trc/ptc/tup/hid/oc/Nrt/ct/Wnt2/Ptp4E/Ptp99A/Ptp10D/Galphaq/Ten-m/mew/pros/zfh1/rho/hh/N/fru/sim/otk/ci/tok/eIF2beta/Mbs/Sox14/robo1/nvy/dac/comm/Wnt5/hep/Cdc42/SrpRalpha/Sdc/trn/Wnt4/tutl/kermit/sqz/TrpRS/msn/Sema1a/babo/fra/Med/Moe/Tig/wts/Hem/beat-Ia/shot/Cdk5/Rac2/Rel/Rho1/bnl/cher/gcm/mirr/Pi3K92E/Shc/chn/drl/tap/CadN/Dref/Nak/NetA/NetB/Rab6/Vang/ara/nompA/side/fz2/gcm2/ko/Rfx/Nrk/Dad/chb/nudC/Tor/Pka-R2/caps/RhoGEF2/wit/WASp/stan/Optix/DAAM/Rok/RhoGAPp190/jim/Cdk5alpha/jbug/bdl/Rich/beat-Ic/Fmr1/seq/nerfin-1/SoxN/ssh/lva/pod1/Dlic/CG11584/hiw/HUWE1/Efhc1.1/CG6867/wgn/rictor/erm/CG14043/BuGZ/Oli/CdGAPr/Dscam1/Mmp2/Drl-2/CG6701/ckn/unc-5/Toll-7/Liprin-gamma/babos/CG13531/spz5/RhoGEF64C/frac/path/Toll-6/sff/wnd/CG8765/Cont/kat-60L1/beat-Vc/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/Dcr-1/atl/XNP/beat-VI/DIP-gamma/RhoGAP100F/Vav/SMC1/Alk/robo3/exex/comm2/dlp/SCAR/ham/Liprin-alpha/CG30456/Unc-115a/Trim9/gogo/LRP1/CG33111/sav/dpr11/Mical/RunxA/nub/spri/sif/Snoo/jing/jeb/vvl/psidin/kra/beta-Spec/alpha-Spec/futsch/Mob4/crb/mmy/cwo/kuz/Mp/mfas/Dscam3/haf/NT1/wb/pdm3/Dbx/Ank2/Dhc64C/Bsg/mid/Ire1/Ptpmeg/CadN2/Sec23/Src64B/cnc/Acsl/Dscam4/dar1/hts/Dark/sick/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/Dscam2/Mtp/tau/Ten-a/unc-104/SelR/Frl/ap/br/S6k/InR/Sac1/Sema5c/trol/btl/sbb/cv-c 276
GO:0120039 plasma membrane bounded cell projection morphogenesis 276/3669 438/8880 0 0 0 acj6/Appl/brm/bsk/chif/chic/eya/D/dnc/EcR/ems/en/E(spl)m8-HLH/Fas2/Fas3/FER/Glt/grn/if/kay/kn/Arpc1/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/Poxm/ras/rut/Scm/Sh/shn/sm/ss/Tm1/trc/ptc/tup/hid/oc/Nrt/ct/Wnt2/Ptp4E/Ptp99A/Ptp10D/Galphaq/Ten-m/mew/pros/zfh1/rho/hh/N/fru/sim/otk/ci/tok/eIF2beta/Mbs/Sox14/robo1/nvy/dac/comm/Wnt5/hep/Cdc42/SrpRalpha/Sdc/trn/Wnt4/tutl/kermit/sqz/TrpRS/msn/Sema1a/babo/fra/Med/Moe/Tig/wts/Hem/beat-Ia/shot/Cdk5/Rac2/Rel/Rho1/bnl/cher/gcm/mirr/Pi3K92E/Shc/chn/drl/tap/CadN/Dref/Nak/NetA/NetB/Rab6/Vang/ara/nompA/side/fz2/gcm2/ko/Rfx/Nrk/Dad/chb/nudC/Tor/Pka-R2/caps/RhoGEF2/wit/WASp/stan/Optix/DAAM/Rok/RhoGAPp190/jim/Cdk5alpha/jbug/bdl/Rich/beat-Ic/Fmr1/seq/nerfin-1/SoxN/ssh/lva/pod1/Dlic/CG11584/hiw/HUWE1/Efhc1.1/CG6867/wgn/rictor/erm/CG14043/BuGZ/Oli/CdGAPr/Dscam1/Mmp2/Drl-2/CG6701/ckn/unc-5/Toll-7/Liprin-gamma/babos/CG13531/spz5/RhoGEF64C/frac/path/Toll-6/sff/wnd/CG8765/Cont/kat-60L1/beat-Vc/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/Dcr-1/atl/XNP/beat-VI/DIP-gamma/RhoGAP100F/Vav/SMC1/Alk/robo3/exex/comm2/dlp/SCAR/ham/Liprin-alpha/CG30456/Unc-115a/Trim9/gogo/LRP1/CG33111/sav/dpr11/Mical/RunxA/nub/spri/sif/Snoo/jing/jeb/vvl/psidin/kra/beta-Spec/alpha-Spec/futsch/Mob4/crb/mmy/cwo/kuz/Mp/mfas/Dscam3/haf/NT1/wb/pdm3/Dbx/Ank2/Dhc64C/Bsg/mid/Ire1/Ptpmeg/CadN2/Sec23/Src64B/cnc/Acsl/Dscam4/dar1/hts/Dark/sick/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/Dscam2/Mtp/tau/Ten-a/unc-104/SelR/Frl/ap/br/S6k/InR/Sac1/Sema5c/trol/btl/sbb/cv-c 276
GO:0000904 cell morphogenesis involved in differentiation 294/3669 474/8880 0 0 0 acj6/Appl/baz/brm/bsk/chif/chic/eya/D/dnc/EcR/ems/en/E(spl)m8-HLH/Fas2/Fas3/FER/Glt/grn/if/kay/kn/Arpc1/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/Poxm/ras/rut/Scm/Sh/shn/sm/ss/Tm1/trc/ptc/tup/hid/oc/Nrt/ct/Wnt2/Ptp4E/Ptp99A/Ptp10D/Galphaq/Ten-m/mew/pros/zfh1/rho/N/fru/sim/otk/ci/ph-p/tok/eIF2beta/Sox14/robo1/nvy/dac/comm/Wnt5/alpha-Cat/hep/Cdc42/SrpRalpha/Sdc/trn/Wnt4/tutl/kermit/sqz/TrpRS/msn/Sema1a/babo/fra/Med/Moe/Tig/wts/Hem/beat-Ia/shot/Cdk5/Rac2/Rel/Rho1/cher/gcm/mirr/Pi3K92E/Shc/chn/drl/tap/CadN/Dref/Nak/NetA/NetB/Rab6/Vang/ara/nompA/side/fz2/gcm2/stck/ko/Rfx/Nrk/Dad/chb/nudC/Tor/Pka-R2/caps/RhoGEF2/wit/WASp/stan/Optix/DAAM/Rok/Traf4/RhoGAPp190/jim/Cdk5alpha/jbug/Ilk/bdl/Rich/beat-Ic/Fmr1/seq/nerfin-1/SoxN/ssh/lva/pigs/pod1/Dlic/CG11584/hiw/HUWE1/Efhc1.1/CG6867/wgn/rictor/erm/CG14043/Megf8/BuGZ/Oli/CdGAPr/Dscam1/Mmp2/Drl-2/CG6701/ckn/unc-5/Toll-7/Liprin-gamma/babos/CG13531/spz5/RhoGEF64C/frac/path/Toll-6/sff/wnd/CG8765/zye/Cont/kat-60L1/beat-Vc/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/Dcr-1/atl/XNP/beat-VI/DIP-gamma/neo/mey/nyo/RhoGAP100F/Vav/SMC1/Alk/robo3/exex/comm2/dlp/SCAR/ham/Liprin-alpha/CG30456/Unc-115a/smog/Trim9/gogo/parvin/LRP1/sav/dpr11/Mical/dyl/Patj/RunxA/Shrm/nub/spri/sif/Snoo/jing/jeb/vvl/psidin/kra/beta-Spec/alpha-Spec/futsch/grh/crb/mmy/cwo/kuz/Mp/mfas/Dscam3/Ppcs/haf/NT1/wb/pdm3/Dbx/Ank2/Dhc64C/LanB1/Bsg/mid/Ire1/Ptpmeg/CadN2/Sec23/Src64B/cnc/Acsl/Dscam4/dar1/hts/Dark/sick/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/Dscam2/Fhos/Mtp/tau/Ten-a/unc-104/LanB2/SelR/Frl/ap/br/S6k/InR/Sac1/Sema5c/trol/tyn/btl/sbb/cv-c 294
GO:0048812 neuron projection morphogenesis 274/3669 436/8880 0 0 0 acj6/Appl/brm/bsk/chif/chic/eya/D/dnc/EcR/ems/en/E(spl)m8-HLH/Fas2/Fas3/FER/Glt/grn/if/kay/kn/Arpc1/Sos/Lim3/sas/LanA/robo2/sens/peb/pk/Poxm/ras/rut/Scm/Sh/shn/sm/ss/Tm1/trc/ptc/tup/hid/oc/Nrt/ct/Wnt2/Ptp4E/Ptp99A/Ptp10D/Galphaq/Ten-m/mew/pros/zfh1/rho/N/fru/sim/otk/ci/tok/eIF2beta/Mbs/Sox14/robo1/nvy/dac/comm/Wnt5/hep/Cdc42/SrpRalpha/Sdc/trn/Wnt4/tutl/kermit/sqz/TrpRS/msn/Sema1a/babo/fra/Med/Moe/Tig/wts/Hem/beat-Ia/shot/Cdk5/Rac2/Rel/Rho1/cher/gcm/mirr/Pi3K92E/Shc/chn/drl/tap/CadN/Dref/Nak/NetA/NetB/Rab6/Vang/ara/nompA/side/fz2/gcm2/ko/Rfx/Nrk/Dad/chb/nudC/Tor/Pka-R2/caps/RhoGEF2/wit/WASp/stan/Optix/DAAM/Rok/RhoGAPp190/jim/Cdk5alpha/jbug/bdl/Rich/beat-Ic/Fmr1/seq/nerfin-1/SoxN/ssh/lva/pod1/Dlic/CG11584/hiw/HUWE1/Efhc1.1/CG6867/wgn/rictor/erm/CG14043/BuGZ/Oli/CdGAPr/Dscam1/Mmp2/Drl-2/CG6701/ckn/unc-5/Toll-7/Liprin-gamma/babos/CG13531/spz5/RhoGEF64C/frac/path/Toll-6/sff/wnd/CG8765/Cont/kat-60L1/beat-Vc/Sra-1/beat-IIb/beat-IIa/RhoGAP93B/Dcr-1/atl/XNP/beat-VI/DIP-gamma/RhoGAP100F/Vav/SMC1/Alk/robo3/exex/comm2/dlp/SCAR/ham/Liprin-alpha/CG30456/Unc-115a/Trim9/gogo/LRP1/CG33111/sav/dpr11/Mical/RunxA/nub/spri/sif/Snoo/jing/jeb/vvl/psidin/kra/beta-Spec/alpha-Spec/futsch/Mob4/crb/mmy/cwo/kuz/Mp/mfas/Dscam3/haf/NT1/wb/pdm3/Dbx/Ank2/Dhc64C/Bsg/mid/Ire1/Ptpmeg/CadN2/Sec23/Src64B/cnc/Acsl/Dscam4/dar1/hts/Dark/sick/dally/dysc/sli/Sema2b/ab/Nrg/CG44153/Dscam2/Mtp/tau/Ten-a/unc-104/SelR/Frl/ap/br/S6k/InR/Sac1/Sema5c/trol/btl/sbb/cv-c 274
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "10 enrichments found"
Code
nrow(up_ekegg)
[1] 10
Code
write.csv(up_ekegg, file = paste0(output_dir, "O1W4KO_vs_ctrl_up_kegg.csv"))
Code
head(up_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme03013 dme03013 Nucleocytoplasmic transport - Drosophila melanogaster (fruit fly) 56/803 99/3488 0.00e+00 0.0000000 0.0000000 mago/RanGAP/lwr/Nup214/Ref1/Mtor/Hel25E/Ran/Nup154/Cbp80/Cbp20/Bin1/Dbp80/Kap-alpha1/Tnpo/mbo/Kap-alpha3/Nup93-1/Nxt1/Aladin/Nup35/Elys/tho2/Tnpo-SR/Nup50/Gle1/tsu/Nup54/Nup62/Rae1/Pym/Upf3/thoc5/thoc7/thoc6/Srrm1/Hpr1/CG7483/RnpS1/Pnn/Ranbp9/Nup93-2/Nup133/Nup98-96/Nup37/Nup358/CG18787/Apl/Ranbp16/Nup153/cdm/Nup160/Nxf3/smt3/ebo/Gp210 56
dme03040 dme03040 Spliceosome - Drosophila melanogaster (fruit fly) 62/803 141/3488 0.00e+00 0.0000009 0.0000008 crn/SmF/l(1)10Bb/mago/snf/Bx42/U2af50/Ref1/Hel25E/noi/Spx/U2af38/ytr/Cbp80/Cbp20/CG2685/l(1)G0007/Prp18/CG4119/x16/hfp/Tango4/CG1622/CG6227/Usp39/tho2/CG3542/Phf5a/CG6686/tsu/CG12343/Cypl/CG2021/CG11586/CG6610/Sf3b6/CG6418/CG10418/Sf3a2/U4-U6-60K/CG14641/Hpr1/CG7483/PQBP1/BCAS2/CG15525/CG17454/c12.1/Prp38/LSm3/Spf45/pea/Prp19/SmD2/SmE/SNRPG/snRNP-U1-C/SmD1/SmB/Cdc5/SF2/caz 62
dme03250 dme03250 Viral life cycle - HIV-1 - Drosophila melanogaster (fruit fly) 20/803 33/3488 3.60e-06 0.0001398 0.0001306 BicD/Fur2/TH1/Snr1/dod/Cdk9/Ran/Gcn5/ear/fal/spt4/Hrs/Tnpo-SR/Cpsf6/Nup358/Jasper/Spt5/Nup153/ALiX/nej 20
dme03082 dme03082 ATP-dependent chromatin remodeling - Drosophila melanogaster (fruit fly) 36/803 81/3488 1.32e-05 0.0003837 0.0003586 Act5C/His2Av/pho/mor/Uch-L5/Snr1/Arp6/HDAC1/Bap60/BCL7-like/CG8677/MRG15/MBD-like/Bap111/CDK2AP1/Arp8/CG31917/Gas41/YL-1/d4/CG10395/MrgBP/CG18004/DMAP1/Eaf6/polybromo/rept/pont/CG12659/Chrac-16/Chrac-14/CG33695/Ino80/Mi-2/Rcd5/Caf1-55 36
dme04120 dme04120 Ubiquitin mediated proteolysis - Drosophila melanogaster (fruit fly) 43/803 107/3488 4.18e-05 0.0008475 0.0007922 ben/fzy/mr/sina/Su(var)2-10/shtd/lwr/eff/Ubc2/Ubc4/Cul1/Cbl/CG2924/SkpA/Cdc16/Ubc10/STUB1/Uba2/Aos1/APC7/UbcE2H/CG4502/CG15237/CG7220/CG7747/Herc4/UbcE2M/Cul6/RhoBTB/Skp2/Fancl/CG2218/ago/Roc2/Ube3a/Ubi-p5E/lmgA/Nedd4/pic/Prp19/Traf6/EloC/Ubc7 43
dme03022 dme03022 Basal transcription factors - Drosophila melanogaster (fruit fly) 21/803 40/3488 4.38e-05 0.0008475 0.0007922 e(y)1/hay/TfIIB/TfIIFalpha/Trf/Taf1/Taf5/TfIIFbeta/TfIIA-L/Taf12/Taf11/Taf2/TfIIA-S/Taf8/Taf7/Taf10b/Tfb4/Taf13/Ssl1/Cdk7/Tfb5 21
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "4 enrichments found"
Code
nrow(down_ekegg)
[1] 4
Code
write.csv(down_ekegg, file = paste0(output_dir, "O1O1W4KO_vs_ctrl_down_kegg.csv"))
Code
head(down_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme04013 dme04013 MAPK signaling pathway - fly - Drosophila melanogaster (fruit fly) 48/1117 93/3488 0.0000562 0.0048255 0.0042684 aop/bi/boss/bsk/chic/dpp/Jra/kay/Sos/peb/Raf/pnt/raw/sev/Btk29A/tkv/Egfr/ttk/vap/vn/hid/RasGAP1/aos/pros/14-3-3zeta/hep/msn/Mef2/phyl/Ask1/Rac2/sty/Shark/Shc/kek1/emb/edl/Mkk4/msk/Tak1/wgn/CYLD/egr/Krn/puc/hkb/Src64B/Duox 48
dme04391 dme04391 Hippo signaling pathway - fly - Drosophila melanogaster (fruit fly) 32/1117 56/3488 0.0000811 0.0048255 0.0042684 baz/bsk/fj/ft/hth/dlg1/sd/tsh/vn/ex/14-3-3zeta/Zyx/wts/jub/wgn/lft/egr/yki/fred/Krn/sav/Patj/crb/app/sdt/d/kibra/Myc/scrib/dally/wg/ds 32
dme00562 dme00562 Inositol phosphate metabolism - Drosophila melanogaster (fruit fly) 25/1117 45/3488 0.0008680 0.0288990 0.0255630 Plc21C/Pi3K59F/Pi3K92E/sktl/Ocrl/CG17896/mtm/Inos/Mipp1/Pten/CG3530/fab1/CG9784/FIG4/PIP5K59B/CG6910/INPP5E/CG17029/Pi4KIIalpha/CG7956/CG42271/5PtaseI/Pi4KIIIalpha/Sac1/IP3K2 25
dme04512 dme04512 ECM-receptor interaction - Drosophila melanogaster (fruit fly) 9/1117 11/3488 0.0009714 0.0288990 0.0255630 Col4a1/LanA/vkg/Hml/CG3168/Tsp/wb/LanB1/LanB2 9

3.8 O2W4_KO

Code
results_O2W4KO_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "O2W4_KO", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_O2W4_KO_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_O2W4KO_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 3772
Code
### MA plot
generate_ma(res)

Code
### Volcano plot
generate_volcano(res, "O2W4_KO", "Ctrl")

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 1770
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "418 enrichments found"
Code
nrow(ego)
[1] 418
Code
write.csv(ego, file = paste0(output_dir, "O2W4KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0051276 chromosome organization 136/1496 298/8880 0 0 0 asp/aub/Cam/CycA/CycB/fs(1)Ya/Hmr/kel/pho/mei-S332/mod(mdg4)/ms(3)K81/mu2/mus301/mus304/Blm/ncd/pim/polo/spn-B/spn-D/spn-E/sub/Su(var)3-7/Su(var)205/Su(var)2-10/thr/Klp61F/Klp67A/shtd/rhi/wapl/piwi/enc/Top1/gammaTub37C/Parp/Irbp/Mtor/prod/glu/Cul1/spel1/HDAC1/pgc/Stat92E/Mcm5/lok/mre11/Mcm7/Nup154/san/SuUR/Bub3/Bap60/Cdc16/Mcm6/cav/Mad1/Nipped-B/Su(var)2-HP2/RecQ5/MRG15/Nup93-1/SMC2/ball/MBD-like/RfC38/APC7/Fign/Pif1/HP6/lid/Nuf2/borr/RfC3/Mcm10/Nup50/Mos/mars/tej/krimp/Nup62/Oxp/Panx/Rae1/rad50/Sse/mad2/nudE/PCID2/Smn/HipHop/Snm1/Hat1/CG8223/BigH1/Mau2/cal1/Fancd2/HP1c/Nup98-96/Mink/cana/RecQ4/cid/Ku80/Chrac-14/SA-2/tefu/Haspin/Rif1/nenya/Nnf1b/His1:CG33801/His1:CG33861/Zwilch/del/tum/egg/Spc25/lmgA/pds5/Su(var)3-3/gwl/RfC4/Incenp/zuc/Psf2/Mi-2/ver/BubR1/Topors/solo/vas/Top2 136
GO:0006259 DNA metabolic process 138/1496 329/8880 0 0 0 aurA/ben/chif/DNaseII/fs(1)Ya/dup/gnu/Hmr/pho/mago/ms(3)K81/mu2/mus101/DNApol-zeta/mus304/Blm/mus312/otu/plu/spn-B/spn-D/spn-E/Su(var)205/Su(var)2-10/Cdk1/tam/Xpc/RpII15/piwi/Gnf1/Top1/PCNA/Parp/CycE/Pms2/RnrL/Irbp/mof/CkIalpha/spel1/tos/smid/Mcm5/lok/mre11/Mcm7/rnh1/Pp4-19C/Orc6/Parg/CG14803/CG14805/CG14814/Drep1/agt/SuUR/Mcm6/cav/Nipped-B/Dek/CG6171/XRCC1/Rad23/RecQ5/Pol31/Ercc1/c(2)M/RfC38/Caf1-180/CG2990/cerv/Arp8/Xrcc2/CG13690/Tfb4/Pif1/HP6/Nse1/Rat1/Mnn1/CG5181/RfC3/Nse4/RPA2/CG9272/Mcm10/CG8679/geminin/mms4/CG9062/Dcr-2/nopo/mip40/DMAP1/rad50/pita/Psf1/Uev1A/Sse/Cdc6/CG3448/ssp/Msh6/HipHop/Snm1/dgrn/MAGE/CG10445/CG9601/Fancl/BigH1/CG5466/Fancd2/RecQ4/mus81/ago/Ku80/Chrac-14/tefu/Rif1/nenya/His1:CG33801/His1:CG33861/AlkB/Ino80/CG42374/pds5/gkt/pic/RfC4/mei-38/Prp19/Psf2/ver/Nup160/Caf1-55/RPA3/Top2 138
GO:0006974 cellular response to DNA damage stimulus 100/1496 233/8880 0 0 0 aurA/ben/mu2/mus101/DNApol-zeta/mus304/Blm/mus312/spn-B/spn-D/Su(var)2-10/Cdk1/tam/Xpc/RpII15/Gnf1/PCNA/Parp/Pms2/Irbp/l(2)dtl/mof/CkIalpha/spel1/tos/Mcm5/lok/Gcn2/mre11/Mcm7/Pp4-19C/Parg/CG14803/CG14814/agt/Mcm6/Nipped-B/Dek/CG6171/XRCC1/Rad23/RecQ5/Ercc1/RfC38/Corp/Caf1-180/cerv/CG8173/Arp8/Xrcc2/CG13690/Tfb4/Pif1/Nse1/Mnn1/CG5181/RfC3/Nse4/RPA2/CG9272/CG8679/mms4/CG9062/DMAP1/rad50/Uev1A/CG3448/Msh6/HipHop/Snm1/CG1218/dgrn/MAGE/CG10445/CG9601/CG11986/Fancl/CG5466/Fancd2/CG7044/RecQ4/mus81/armi/Ku80/Chrac-14/tefu/Claspin/Ino80/CG42374/pds5/gkt/pic/RfC4/Prp19/grp/hpo/Psf2/Nup160/CG43295/RPA3 100
GO:0000278 mitotic cell cycle 160/1496 465/8880 0 0 0 bru1/asp/aub/aurA/ben/BicC/csw/CycA/CycB/dsh/png/fs(1)Ya/dup/fzy/gnu/Hmr/sti/Rcc1/Map205/twe/mor/mud/mus101/mus304/ncd/pim/plu/polo/pum/smo/spn-E/sub/Su(var)205/thr/Cdk1/Klp61F/Klp67A/shtd/wapl/Su(H)/RpII15/PCNA/PpD3/gammaTub37C/CycD/CycJ/CycE/Mad/tsr/Wee1/Arp53D/l(2)dtl/Mtor/prod/noi/glu/CycB3/lok/Ran/mre11/asun/Pp4-19C/Orc6/san/Dlc90F/Bub3/Cdc16/Mcm6/Mad1/Nipped-B/Su(var)2-HP2/Grip75/Grip128/Apc2/RecQ5/SMC2/ball/Pol31/CHORD/Aos1/APC7/Nek2/Aladin/PPP4R2r/CG8173/Rab35/Fign/Mnn1/Nuf2/CG5181/borr/Grip71/Mcm10/geminin/mms4/dgt5/mars/nopo/dgt3/Rae1/JMJD5/Psf1/Sse/mad2/Cdc6/nudE/Syx13/Skp2/PEK/Mau2/cal1/Fancd2/CG13599/Dis3/Mink/CG5003/dgt6/pont/pins/cana/RecQ4/MED22/mus81/cid/ago/tefu/Haspin/Rif1/Nnf1b/Tap42/Claspin/Zwilch/del/ALiX/tum/Spc25/alphaTub67C/lmgA/zld/bora/pds5/gwl/wisp/mei-38/Incenp/grp/nej/aPKC/CG43295/Cdk7/BubR1/Caf1-55/Det/smt3/Smr/Ptp61F/vas/Tlk/Top2/betaTub56D 160
GO:0000280 nuclear division 116/1496 296/8880 0 0 0 bru1/aub/aurA/Axs/cup/CycA/CycB/fzy/Hmr/ifc/sti/twe/mei-S332/mod(mdg4)/mud/mus301/mus304/Blm/mus312/ncd/pim/polo/RanGAP/spn-B/spn-D/spn-E/sub/thr/Cdk1/Klp61F/Klp67A/shtd/Xpc/yem/gammaTub37C/CycE/mtrm/dhd/Pp2B-14D/Mtor/prod/glu/Cul1/mael/Mcm5/Nup154/san/fws/Bub3/Cdc16/Mad1/Nipped-B/Su(var)2-HP2/Grip75/Grip128/RecQ5/SMC2/ball/Ercc1/c(2)M/CHORD/APC7/Nek2/Aladin/Fign/HP6/Nuf2/borr/Grip71/CG8679/mms4/dgt5/Mos/dgt3/Rae1/thoc5/Sse/mad2/Pex2/nudE/Syx13/Mau2/cal1/Fancd2/CG13599/Mink/dgt6/pont/cana/mus81/cid/ago/Haspin/nenya/Nnf1b/Zwilch/tef/tum/sra/Spc25/alphaTub67C/lmgA/bora/pds5/gwl/wisp/mei-38/Incenp/nej/BubR1/Det/mamo/Topors/solo/vas/Top2 116
GO:0048285 organelle fission 119/1496 312/8880 0 0 0 bru1/aub/aurA/Axs/cup/CycA/CycB/fzy/Hmr/ifc/sti/twe/mei-S332/mod(mdg4)/mud/mus301/mus304/Blm/mus312/ncd/pim/polo/RanGAP/spn-B/spn-D/spn-E/sub/thr/Cdk1/Klp61F/Klp67A/shtd/Xpc/yem/gammaTub37C/CycE/mtrm/dhd/Pp2B-14D/Mtor/prod/glu/Cul1/mael/Mcm5/Nup154/Pgam5/san/fws/Bub3/Cdc16/Mad1/Nipped-B/Su(var)2-HP2/Grip75/Grip128/RecQ5/SMC2/ball/Ercc1/c(2)M/CHORD/APC7/Pink1/Nek2/CG10761/Aladin/Fign/HP6/Nuf2/borr/Grip71/CG8679/mms4/dgt5/Mos/dgt3/Rae1/thoc5/Sse/mad2/Pex2/nudE/Syx13/Mau2/cal1/Fancd2/CG13599/Mink/dgt6/pont/cana/mus81/cid/ago/Haspin/nenya/Nnf1b/Zwilch/tef/tum/sra/Spc25/alphaTub67C/lmgA/bora/pds5/gwl/wisp/mei-38/Incenp/nej/BubR1/Det/mamo/Topors/solo/vas/Top2 119
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 4322
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "217 enrichments found"
Code
nrow(ego)
[1] 217
Code
write.csv(ego, file = paste0(output_dir, "O2W4KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0034330 cell junction organization 209/3508 371/8880 0 0e+00 0e+00 Act57B/Appl/arm/awd/baz/Lar/Fas2/Galphai/Galphas/kay/dlg1/Gli/LanA/mr/Atpalpha/put/rut/sax/shi/slo/stau/kst/Ten-m/Fur1/Gad1/Syn/pros/Rap1/gig/ewg/Rab3/comm/brat/Rac1/Cdc42/htl/Sdc/cora/sinu/jar/Nos/repo/nmo/pck/Cdk5/Nrx-IV/Rac2/Rho1/drl/Rab11/RtGEF/lbm/fz2/klg/mGluR/Pdk1/Dad/caps/wit/sns/WASp/trio/Tsp2A/stan/Syx4/moody/SkpA/Srpk79D/unc-13/TBPH/wrapper/Top3beta/Tip60/Rok/drpr/cindr/Ulp1/kirre/Rich/Vap33/Grip/raptor/pod1/dpr14/hiw/Lrp4/Frq1/shv/bark/DIP-eta/Nlg2/dpr19/PICK1/spict/Cul2/kune/Mmp2/Dyb/Syn2/dpr13/Magi/Liprin-gamma/CG13506/PIP5K59B/CG10339/Mmp1/CG13704/crim/RhoGAP68F/Tsf2/CG9628/sff/Ssk/cmpy/M6/CG7166/Syn1/Ada2b/beag/ATP6AP2/dpr5/dpr15/Octbeta2R/dpr9/Sra-1/Lrrk/Nrx-1/Octbeta1R/spas/BOD1/DIP-gamma/Cul5/RhoGAP100F/mav/Cals/Vav/Kat60/Aplip1/dpr1/Prosap/dpr6/Atx2/dlp/cpx/bchs/metro/mesh/Nlg1/CG31475/DIP-theta/DIP-zeta/DIP-kappa/dpr10/dpr8/Synd/dpr11/Mical/dpr4/dpr3/Frq2/Nlg3/Nlg4/shakB/dpr12/DIP-delta/jeb/beta-Spec/futsch/DIP-beta/brp/Mob2/DIP-epsilon/Ca-beta/mtg/Atg1/dpr21/wol/stj/bou/CG42673/Ank2/CG42748/sdt/Tl/Rbp/Rab8/Arp3/Src64B/cac/scrib/nwk/dysc/Fife/CaMKII/Nrg/l(3)neo38/CG44325/Neto/CrebB/ringer/Nsf2/GckIII/stai/kis/btsz/Ten-a/Sac1/wg/Rab1 209
GO:0040011 locomotion 234/3508 429/8880 0 1e-07 1e-07 brm/bsk/chic/Dab/Lar/E(Pc)/Fas2/Fas3/FER/Glt/grk/grn/if/kay/egh/dlg1/Gli/Ca-alpha1D/Catsup/Lim3/sas/LanA/sens/Mlc2/peb/pk/put/sm/spz/svp/Egfr/ttk/tup/vir/vn/w/uzip/Pp1-87B/ct/Ptp99A/Ten-m/mew/pros/hh/fru/Hrb27C/otk/eIF2beta/ple/robo1/nvy/slbo/Orc3/dac/Con/comm/Tbh/Rac1/Cdc42/htl/Nmdar1/Sdc/Wnt4/tutl/kermit/msn/Dop1R1/Moe/Tig/nmo/beat-Ia/shot/Cdk5/Gycalpha99B/Ptp69D/Rac2/Rho1/Tie/bnl/cher/mirr/tim/drl/tap/CadN/NetA/NetB/Rab6/Vang/fz2/14-3-3epsilon/ko/dom/Nrk/Dad/CG7179/chb/Pka-R2/caps/Pdf/CG3857/trio/stan/DAAM/PlexA/homer/TBPH/Timp/Stam/Cdk5alpha/jbug/Rich/nerfin-1/raptor/pod1/Clic/hiw/CG9784/CG6867/wgn/Nup205/daw/Pvf2/Pvr/Bug22/CG17118/Oli/Dscam1/Mmp2/ths/Drl-2/ckn/Hr51/unc-5/Syn2/CG6805/Toll-7/Lrt/spz5/RhoGEF64C/atk/M6/Syn1/laza/ATP6AP2/beat-Vc/beat-Va/beat-Vb/Sra-1/beat-IIb/beat-IIa/Lrrk/Ktl/RhoGAP93B/Nrx-1/AP-1sigma/XNP/CG6330/beat-VI/DIP-gamma/RhoGAP100F/Vav/Kat60/boi/tai/robo3/dlp/SCAR/Proc/CG30456/Unc-115a/Trim9/Krn/gogo/shep/dpr11/Mical/pain/CG34353/Pvf3/spri/jeb/vvl/l(3)L1231/sls/puc/kra/beta-Spec/meigo/mmy/kuz/flr/Diap1/Mp/Pdfr/Dscam3/conv/haf/NT1/wb/Msp300/mid/Ptpmeg/CadN2/Tl/Src64B/Evi5/Acsl/hts/Hmgcr/dally/Stacl/dysc/sli/Sema2b/ab/Mhc/Nrg/CG44153/Neto/beat-IIIa/rdgC/ttv/lov/kis/Sec15/Ten-a/ap/br/Sac1/Sema5c/trol/sbb/Rab1 234
GO:0099537 trans-synaptic signaling 168/3508 292/8880 0 1e-07 1e-07 Ace/nAChRalpha1/mAChR-A/nAChRbeta1/nAChRalpha2/boss/Cat/Pka-C1/ChAT/comt/dnc/Fas2/Galphas/dlg1/mr/na/Atpalpha/rt/rut/Sh/shi/slo/svp/nAChRbeta2/5-HT1A/Csp/Syt1/Rdl/Ten-m/5-HT7/Rop/Syn/Cbp53E/Plc21C/TkR99D/Rap1/Shal/Rab3/CanB/Lcch3/Cdc42/Nmdar1/Dcp-1/Snap25/Dop1R1/nSyb/Syx1A/Arf51F/Ipp/stnB/mGluR/dtr/wit/b6/Oamb/GluClalpha/Syx4/Srpk79D/unc-13/TBPH/Tip60/X11L/Rep/Rab3-GAP/GABA-B-R2/GluRIID/nAChRalpha5/NtR/Vap33/hiw/Rab3-GEF/Frq1/wgn/VGlut/Snapin/Nlg2/CG7781/CG14275/nAChRalpha6/atilla/twit/CG9336/CG9338/Dyb/CG6329/Syngr/VGAT/DAT/P32/CG13531/CG10339/spz5/FMRFaR/unc-13-4A/PGRP-LC/crim/rogdi/SecCl/SLC22A/Neurochondrin/beag/ohgt/DCAF12/KaiR1D/Nrx-1/spas/CDase/kek6/Syt7/Ekar/CG11155/Cals/cpx/GluRIIC/RIC-3/Nlg1/GluRIIE/CG31475/CG31675/CG31676/CG31935/Csas/CG32447/X11Lbeta/Rim/Cadps/Frq2/Nlg3/Nlg4/jeb/beta-Spec/alphaSnap/brp/Ca-beta/qvr/Atg1/stj/Syt12/Sytbeta/rtv/NT1/Exn/nemy/Ank2/Rbp/Shab/Src64B/SERCA/cac/5-HT1B/Acsl/Ih/nwk/cpo/dysc/Slob/Fife/CaMKII/Neto/nAChRalpha4/Nsf2/Sec5/Sec6/Snap24/Esyt2/Ten-a/VAChT/cv-c 168
GO:0003008 system process 252/3508 472/8880 0 1e-07 1e-07 nAChRalpha1/nAChRbeta1/nAChRalpha2/Appl/arm/boss/bsk/bw/Pka-C1/salr/Col4a1/dnc/eag/EcR/Fas2/Galphas/Hn/if/Ca-alpha1D/Mlc1/Mlc2/na/Atpalpha/Pkc98E/pnr/Rh4/rut/Sh/shi/Btk29A/ss/stau/sws/Egfr/vn/w/Pp1-87B/Tm2/nAChRbeta2/5-HT1A/up/ct/Rdl/Klp64D/mew/Mdr65/Oct-TyrR/Gad1/Syn/pros/TkR99D/svr/hop/14-3-3zeta/Shal/Orc3/Galphaf/Lcch3/Tbh/Rac1/Cdc42/Nmdar1/cora/sinu/dia/Obp28a/RyR/Ms/Dop1R1/ine/Mvl/Nos/Fdh/Nrx-IV/cher/tim/drl/Drip/nompA/vri/H15/Dyrk2/klg/mGluR/ppk/geko/Rfx/CRMP/Clk/teq/CG3078/sns/b6/Oamb/GluClalpha/Scgdelta/Amnionless/Eb1/Inx7/Amph/CG4928/kirre/Kap3/Rh50/nAChRalpha5/aru/NtR/Fkbp59/CG5921/dpr14/fh/Uba5/Sclp/hiw/Mec2/CG8028/Ir25a/meng/Nha1/PrBP/nAChRalpha6/dpr19/AstC/MRP/nrv3/Dscam1/Obp44a/KCNQ/CAP/ClC-b/Sans/dpr13/Oatp58Dc/Eglp2/Eglp4/ITP/mri/Chd64/DopEcR/ppk26/Bmcp/INPP5E/Or71a/SecCl/SLC22A/CG14636/NPFR/Neurochondrin/ATP6AP2/mura/dpr5/ATP8B/dpr15/dpr9/Sdhaf3/Nrx-1/CG6972/Root/sosie/PQBP1/BOD1/Obp99a/Obp99c/Obp99d/Obp99b/Cad99C/CG15537/G9a/dpr1/dpr6/robo3/lilli/Atx2/HDAC4/prt/Proc/Obp83g/Obp83ef/Obp83cd/SKIP/Gba1b/wry/CG31935/dpr10/Rbfox1/Krn/dpr8/Cubn/dpr11/dpr4/Nmdar2/dpr3/Dop2R/Unc-89/pain/Ctr1A/Or19b/dbr/dpr12/Nox/NaCP60E/sif/kra/DCX-EMAP/futsch/brp/Mob2/drd/qvr/Mp/2mit/Atg1/stj/salm/nemy/Npc1b/Ank2/dikar/mid/f/Shab/Src64B/AGO1/Ptr/SERCA/cac/Fic/scrib/Nha2/spoon/dysc/orb2/CG43778/Ca-alpha1T/CaMKII/Mhc/Rgk1/tn/rad/CrebB/rdgC/nAChRalpha4/SclA/SclB/kis/ebo/Ggamma30A/Tmc/wupA/Phm/wg/para 252
GO:0007268 chemical synaptic transmission 167/3508 291/8880 0 1e-07 1e-07 Ace/nAChRalpha1/mAChR-A/nAChRbeta1/nAChRalpha2/boss/Cat/Pka-C1/ChAT/comt/dnc/Fas2/Galphas/dlg1/mr/na/Atpalpha/rt/rut/Sh/shi/slo/svp/nAChRbeta2/5-HT1A/Csp/Syt1/Rdl/Ten-m/5-HT7/Rop/Syn/Cbp53E/Plc21C/TkR99D/Rap1/Shal/Rab3/CanB/Lcch3/Cdc42/Nmdar1/Dcp-1/Snap25/Dop1R1/nSyb/Syx1A/Arf51F/Ipp/stnB/mGluR/dtr/wit/b6/Oamb/GluClalpha/Syx4/Srpk79D/unc-13/TBPH/Tip60/X11L/Rep/Rab3-GAP/GABA-B-R2/GluRIID/nAChRalpha5/NtR/Vap33/hiw/Rab3-GEF/Frq1/wgn/VGlut/Snapin/Nlg2/CG7781/CG14275/nAChRalpha6/atilla/twit/CG9336/CG9338/Dyb/CG6329/Syngr/VGAT/DAT/P32/CG13531/CG10339/spz5/FMRFaR/unc-13-4A/PGRP-LC/crim/rogdi/SecCl/SLC22A/Neurochondrin/beag/ohgt/DCAF12/KaiR1D/Nrx-1/spas/CDase/kek6/Syt7/Ekar/CG11155/Cals/cpx/GluRIIC/RIC-3/Nlg1/GluRIIE/CG31675/CG31676/CG31935/Csas/CG32447/X11Lbeta/Rim/Cadps/Frq2/Nlg3/Nlg4/jeb/beta-Spec/alphaSnap/brp/Ca-beta/qvr/Atg1/stj/Syt12/Sytbeta/rtv/NT1/Exn/nemy/Ank2/Rbp/Shab/Src64B/SERCA/cac/5-HT1B/Acsl/Ih/nwk/cpo/dysc/Slob/Fife/CaMKII/Neto/nAChRalpha4/Nsf2/Sec5/Sec6/Snap24/Esyt2/Ten-a/VAChT/cv-c 167
GO:0098916 anterograde trans-synaptic signaling 167/3508 291/8880 0 1e-07 1e-07 Ace/nAChRalpha1/mAChR-A/nAChRbeta1/nAChRalpha2/boss/Cat/Pka-C1/ChAT/comt/dnc/Fas2/Galphas/dlg1/mr/na/Atpalpha/rt/rut/Sh/shi/slo/svp/nAChRbeta2/5-HT1A/Csp/Syt1/Rdl/Ten-m/5-HT7/Rop/Syn/Cbp53E/Plc21C/TkR99D/Rap1/Shal/Rab3/CanB/Lcch3/Cdc42/Nmdar1/Dcp-1/Snap25/Dop1R1/nSyb/Syx1A/Arf51F/Ipp/stnB/mGluR/dtr/wit/b6/Oamb/GluClalpha/Syx4/Srpk79D/unc-13/TBPH/Tip60/X11L/Rep/Rab3-GAP/GABA-B-R2/GluRIID/nAChRalpha5/NtR/Vap33/hiw/Rab3-GEF/Frq1/wgn/VGlut/Snapin/Nlg2/CG7781/CG14275/nAChRalpha6/atilla/twit/CG9336/CG9338/Dyb/CG6329/Syngr/VGAT/DAT/P32/CG13531/CG10339/spz5/FMRFaR/unc-13-4A/PGRP-LC/crim/rogdi/SecCl/SLC22A/Neurochondrin/beag/ohgt/DCAF12/KaiR1D/Nrx-1/spas/CDase/kek6/Syt7/Ekar/CG11155/Cals/cpx/GluRIIC/RIC-3/Nlg1/GluRIIE/CG31675/CG31676/CG31935/Csas/CG32447/X11Lbeta/Rim/Cadps/Frq2/Nlg3/Nlg4/jeb/beta-Spec/alphaSnap/brp/Ca-beta/qvr/Atg1/stj/Syt12/Sytbeta/rtv/NT1/Exn/nemy/Ank2/Rbp/Shab/Src64B/SERCA/cac/5-HT1B/Acsl/Ih/nwk/cpo/dysc/Slob/Fife/CaMKII/Neto/nAChRalpha4/Nsf2/Sec5/Sec6/Snap24/Esyt2/Ten-a/VAChT/cv-c 167
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "17 enrichments found"
Code
nrow(up_ekegg)
[1] 17
Code
write.csv(up_ekegg, file = paste0(output_dir, "O2W4KO_vs_ctrl_up_kegg.csv"))
Code
head(up_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme03013 dme03013 Nucleocytoplasmic transport - Drosophila melanogaster (fruit fly) 52/548 99/3488 0.0000000 0.0000000 0.0000000 mago/RanGAP/Nup214/Ref1/Mtor/Hel25E/Ran/emb/Nup154/Cse1/Cbp80/Cbp20/Bin1/Dbp80/Kap-alpha1/Tnpo/mbo/msk/Kap-alpha3/Nup93-1/Nxt1/Aladin/Nup35/Elys/tho2/Tnpo-SR/Nup50/Gle1/Nup54/Nup62/Rae1/Pym/thoc5/thoc7/Srrm1/Hpr1/Pnn/Ranbp9/Nup98-96/Nup37/Nup358/CG18787/Arts/Apl/Ranbp16/Nup153/Snup/cdm/Nup160/Nxf3/smt3/Gp210 52
dme04330 dme04330 Notch signaling pathway - Drosophila melanogaster (fruit fly) 14/548 25/3488 0.0000040 0.0002174 0.0001872 dsh/gro/H/Su(H)/Bx42/HDAC1/Gcn5/Atx-1/aph-1/CG6843/Nct/pen-2/nej/Psn 14
dme03250 dme03250 Viral life cycle - HIV-1 - Drosophila melanogaster (fruit fly) 16/548 33/3488 0.0000095 0.0003495 0.0003010 BicD/Fur2/Snr1/Cdk9/Ran/Gcn5/emb/ear/fal/Hrs/Tnpo-SR/Cpsf6/Nup358/Nup153/ALiX/nej 16
dme04120 dme04120 Ubiquitin mediated proteolysis - Drosophila melanogaster (fruit fly) 34/548 107/3488 0.0000193 0.0005309 0.0004572 ben/fzy/sina/Su(var)2-10/shtd/Ubc2/Ubc4/Cul1/Cbl/CG2924/Cdc16/Ubc10/STUB1/Uba2/Aos1/APC7/UbcE2H/CG4502/CG7220/Cul6/RhoBTB/Skp2/Fancl/CG2218/ago/Roc2/Ube3a/Ubi-p5E/lmgA/Nedd4/pic/Prp19/Traf6/EloC 34
dme03410 dme03410 Base excision repair - Drosophila melanogaster (fruit fly) 14/548 30/3488 0.0000598 0.0013161 0.0011335 tam/Gnf1/PCNA/Parp/Parg/XRCC1/Pol31/RfC38/RfC3/CG9272/CG9601/Chrac-14/gkt/RfC4 14
dme03082 dme03082 ATP-dependent chromatin remodeling - Drosophila melanogaster (fruit fly) 26/548 81/3488 0.0001546 0.0028341 0.0024408 His2Av/pho/mor/Uch-L5/Iswi/Snr1/HDAC1/Bap60/CG8677/MRG15/MBD-like/Bap111/Arp8/d4/CG10395/CG18004/DMAP1/polybromo/rept/pont/CG12659/Chrac-14/CG33695/Ino80/Mi-2/Caf1-55 26
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "2 enrichments found"
Code
nrow(down_ekegg)
[1] 2
Code
write.csv(down_ekegg, file = paste0(output_dir, "O2W4KO_vs_ctrl_down_kegg.csv"))
Code
head(down_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme04080 dme04080 Neuroactive ligand-receptor interaction - Drosophila melanogaster (fruit fly) 35/1072 67/3488 0.000173 0.0205865 0.0187564 mAChR-A/baz/alphaTry/5-HT1A/Rdl/5-HT7/GluRIA/TkR99D/Lcch3/betaTry/deltaTry/gammaTry/Nmdar1/epsilonTry/Dop1R1/Lgr1/mGluR/AkhR/GABA-B-R2/CCKLR-17D3/GABA-B-R3/CCHa2-R/Dh44-R2/Dh44-R1/CrzR/CCAP-R/kappaTry/CCHa1-R/Dh31-R/Nmdar2/Dop2R/Gpb5/GABA-B-R1/5-HT2B/5-HT1B 35
dme04391 dme04391 Hippo signaling pathway - fly - Drosophila melanogaster (fruit fly) 29/1072 56/3488 0.000746 0.0443890 0.0404429 Act5C/baz/bsk/ft/hth/dlg1/sd/vn/mts/ex/tws/14-3-3zeta/dia/wts/14-3-3epsilon/wgn/Dlg5/mats/p53/fred/Krn/app/sdt/d/Myc/scrib/dally/wg/ds 29

3.9 O1O2_KO

Code
results_O1O2KO_CTRL = get_dds_res(wald_dds, 
                             contrast = c("condition", "O1O2_KO", "Ctrl"),
                             ensembl.genes, 
                             shrink = TRUE)
[1] "condition_O1O2_KO_vs_Ctrl"
using 'apeglm' for LFC shrinkage. If used in published research, please cite:
    Zhu, A., Ibrahim, J.G., Love, M.I. (2018) Heavy-tailed prior distributions for
    sequence count data: removing the noise and preserving large differences.
    Bioinformatics. https://doi.org/10.1093/bioinformatics/bty895
Code
res = results_O1O2KO_CTRL 

hist(res$pvalue, xlab = "p-value", ylab = "Frequency")

Number of significant genes (padj < 0.1 & log2FoldChange < log2(1.5)):

Code
sum(res$padj < 0.1 & abs(res$log2FoldChange) > lfc.threshold, na.rm = T)
[1] 3124
Code
### MA plot
generate_ma(res)

Code
### Volcano plot
generate_volcano(res, "O1O2_KO", "Ctrl")

Code
df_res = as.data.frame(res)
# Drop rows with NAs in pvalue - these are genes with high variability
df_res = df_res[complete.cases(df_res[ ,c("pvalue", "padj")]), ] 
Code
up_genes <- rownames(res[res$padj < 0.1 &
                      res$log2FoldChange > lfc.threshold, ])

length(up_genes)
[1] 1183
Code
ego <- plotEGO_dm(up_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, upregulated genes")
Running GO for organism = drosophila melanogaster

[1] "194 enrichments found"
Code
nrow(ego)
[1] 194
Code
write.csv(ego, file = paste0(output_dir, "O1O2KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0006259 DNA metabolic process 97/976 329/8880 0 0 0 aurA/fs(1)Ya/gnu/hay/Hmr/mei-9/mago/mu2/mus101/mus304/otu/rec/spn-B/spn-D/spn-E/Su(var)205/Cdk1/Xpac/piwi/Gnf1/PCNA/DNApol-alpha73/CycE/Pms2/RnrL/Irbp/mof/spel1/tos/Mcm5/lok/mre11/Orc6/CG14803/CG14814/Drep1/agt/Rad17/Cdc45/CG6171/Rad1/CG10254/Pol31/Ercc1/c(2)M/CG4078/CG2990/Psf3/cerv/CG8142/Arp8/Xrcc2/CG13690/Pif1/HP6/Nse1/RfC3/RPA2/CG9272/CG8679/mms4/CG9062/Ctf4/nopo/rad50/Psf1/Uev1A/Sse/CG3448/Msh6/HipHop/Snm1/dgrn/CG10445/CG9601/BigH1/CG6985/Sld5/RecQ4/mus81/Ku80/Chrac-14/tefu/pch2/Tengl3/CG31807/hd/gkt/RfC4/mrn/Xpd/Psf2/ver/Caf1-55/RPA3/Pol32/vilya 97
GO:0006974 cellular response to DNA damage stimulus 77/976 233/8880 0 0 0 aurA/hay/mei-9/mu2/mus101/mus304/spn-B/spn-D/Cdk1/Xpac/Gnf1/PCNA/Pms2/Irbp/l(2)dtl/mof/spel1/tos/Mcm5/lok/mre11/CG14803/CG14814/agt/Rad17/Cdc45/CG6171/Rad1/CG10254/Ercc1/CG4078/Corp/cerv/CG8173/CG8142/Arp8/Xrcc2/CG9578/CG13690/Pif1/Nse1/RfC3/RPA2/CG9272/CG8679/mms4/CG9062/Ctf4/rad50/Uev1A/CG3448/Msh6/HipHop/Snm1/CG1218/dgrn/CG10445/CG9601/CG11986/CG6985/Sld5/RecQ4/mus81/armi/Ku80/Chrac-14/tefu/CG31807/Claspin/gkt/RfC4/mrn/grp/Xpd/Psf2/RPA3/Pol32 77
GO:0006281 DNA repair 64/976 186/8880 0 0 0 aurA/hay/mei-9/mu2/mus304/spn-B/spn-D/Xpac/Gnf1/PCNA/Pms2/Irbp/mof/spel1/tos/Mcm5/lok/mre11/CG14803/CG14814/agt/Rad17/Cdc45/CG6171/Rad1/CG10254/Ercc1/CG4078/cerv/CG8142/Arp8/Xrcc2/CG13690/Pif1/Nse1/RfC3/RPA2/CG9272/CG8679/mms4/CG9062/Ctf4/rad50/Uev1A/CG3448/Msh6/Snm1/dgrn/CG10445/CG9601/CG6985/Sld5/RecQ4/mus81/Ku80/tefu/CG31807/gkt/RfC4/mrn/Xpd/Psf2/RPA3/Pol32 64
GO:0051276 chromosome organization 78/976 298/8880 0 0 0 aub/CycA/CycB/fs(1)Ya/Hmr/mu2/mus301/mus304/rec/spn-B/spn-D/spn-E/sub/Su(var)3-7/Su(var)205/Klp67A/shtd/rhi/wapl/Xpac/piwi/gammaTub37C/Irbp/prod/spel1/pgc/Mcm5/lok/mre11/san/Cdc45/MRG15/ball/CG4078/APC7/Fign/Pif1/HP6/lid/RfC3/Nup50/Mos/tej/krimp/Panx/Rae1/rad50/Sse/mad2/PCID2/Smn/HipHop/Snm1/Hat1/CG8223/BigH1/cal1/Nup133/HP1c/RecQ4/cid/Ku80/Chrac-14/SA-2/tefu/Haspin/Nnf1b/Zwilch/del/gwl/RfC4/Xpd/Psf2/ver/BubR1/solo/vas/vilya 78
GO:0006310 DNA recombination 40/976 106/8880 0 0 0 mei-9/mago/mus304/rec/spn-B/spn-D/spn-E/piwi/spel1/tos/Mcm5/mre11/CG14803/Cdc45/Ercc1/c(2)M/CG4078/cerv/Xrcc2/Pif1/Nse1/RPA2/CG8679/mms4/CG9062/rad50/CG3448/dgrn/BigH1/CG6985/Sld5/RecQ4/mus81/Ku80/pch2/Xpd/Psf2/RPA3/Pol32/vilya 40
GO:0051321 meiotic cell cycle 61/976 215/8880 0 0 0 aub/Axs/cort/cup/twe/mei-9/mus301/mus304/rec/spn-B/spn-D/spn-E/sub/Cdk1/Klp67A/yem/gammaTub37C/mtrm/dhd/Pp2B-14D/mael/Mcm5/lok/Ran/mre11/mei-P26/Grip84/Grip75/Grip128/ball/Ercc1/c(2)M/HP6/lid/Grip71/CG8679/mms4/Mos/tej/krimp/Rae1/rad50/wac/Sse/Syx13/Msh6/cal1/CG6985/mus81/tefu/pch2/del/tef/sra/gwl/wisp/BubR1/solo/vas/vilya/elgi 61
Code
down_genes <- rownames(df_res[df_res$padj < 0.1 &
                      df_res$log2FoldChange < lfc.threshold, ])

length(down_genes)
[1] 4119
Code
ego <- plotEGO_dm(down_genes, universe = rownames(dds), 
               ont = "BP",
               title = "GO, downregulated genes")
Running GO for organism = drosophila melanogaster

[1] "127 enrichments found"
Code
nrow(ego)
[1] 127
Code
write.csv(ego, file = paste0(output_dir, "O1O2KO_vs_ctrl_down_go.csv"))
Code
head(ego)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
GO:0040003 chitin-based cuticle development 105/3316 147/8880 0e+00 0.00e+00 0.0e+00 amd/arm/ck/ecd/Gld/kkv/knk/l(2)gl/ttk/CrebA/Ccp84Ag/Ccp84Ae/Ccp84Ad/Ccp84Ac/Ccp84Ab/Ccp84Aa/jar/Acp1/fz2/Idgf1/pio/Lcp65Ag2/Lcp65Ag2/Lcp65Af/Lcp65Ac/Lcp65Ab2/Alas/Gasp/Axn/Vajk3/Cpr35B/Vajk1/SoxN/TwdlT/Cpr11A/obst-A/Cpr23B/obst-E/TwdlE/Vajk2/Hakai/Twdlbeta/Cpr49Ac/Cpr49Ae/Cpr49Ag/Cpr50Cb/Cpr51A/resilin/Cpr56F/Cpr62Bb/Cpr62Bc/Cht7/Cpr64Ab/Cpr64Ad/Cpr65Ea/Cpr65Eb/Cpr65Ec/Cpr66Cb/Cpr67Fa1/Cpr67Fa2/Cpr76Bc/Cpr78Cc/TwdlF/TwdlG/TwdlV/Cht5/TwdlW/Cpr92F/TwdlM/TwdlP/TwdlL/TwdlN/TwdlD/TwdlQ/TwdlC/Cpr97Ea/Cpr97Eb/Cpr100A/Cpr65Ax2/Vajk4/TwdlH/fred/Cpr66D/Cpr65Av/TwdlY/TwdlX/dpy/Cpr31A/Idgf5/pan/Cpr47Eg/puc/Tb/CG34461/pot/mmy/drd/serp/wol/verm/dsx-c73A/sdt/Sec24CD/wg/stw 105
GO:0042335 cuticle development 122/3316 185/8880 0e+00 0.00e+00 0.0e+00 amd/arm/ck/Ddc/ect/ecd/Gld/kkv/knk/l(2)gl/sha/ttk/Ser/CrebA/Ccp84Ag/Ccp84Ae/Ccp84Ad/Ccp84Ac/Ccp84Ab/Ccp84Aa/ple/Cyp4g1/jar/fng/Acp1/p38a/fz2/Idgf1/pio/Lcp65Ag2/Lcp65Ag2/Lcp65Af/Lcp65Ac/Lcp65Ab2/Alas/WASp/Mkk4/Gasp/gus/Axn/Vajk3/Cpr35B/Vajk1/SoxN/TwdlT/Cpr11A/obst-A/Cpr23B/snsl/obst-E/TwdlE/Vajk2/Hakai/Twdlbeta/Cpr49Ac/Cpr49Ae/Cpr49Ag/Cpr50Cb/Cpr51A/resilin/Cpr56F/Cpr62Bb/Cpr62Bc/Cht7/Cpr64Ab/Cpr64Ad/Cpr65Ea/Cpr65Eb/Cpr65Ec/Cpr66Cb/Cpr67Fa1/Cpr67Fa2/Cpr76Bc/Cpr78Cc/TwdlF/TwdlG/TwdlV/Cht5/TwdlW/Cpr92F/TwdlM/TwdlP/TwdlL/TwdlN/TwdlD/TwdlQ/TwdlC/Cpr97Ea/Cpr97Eb/snu/Cpr100A/boi/Cpr65Ax2/FASN2/Vajk4/TwdlH/fred/Cpr66D/Cpr65Av/TwdlY/TwdlX/dpy/Cpr31A/Idgf5/pan/Cpr47Eg/Desat1/puc/Tb/CG34461/pot/mmy/drd/serp/wol/verm/dsx-c73A/sdt/f/Sec24CD/wg/stw 122
GO:0003008 system process 241/3316 472/8880 0e+00 3.00e-07 3.0e-07 nAChRalpha1/nAChRbeta1/nAChRalpha2/Appl/arm/boss/bsk/bw/Cam/salr/Col4a1/ck/Ddc/dnc/eag/for/Galphas/Hn/if/inaD/Ca-alpha1D/dco/Mlc1/Mlc2/na/Atpalpha/Pkc98E/pnr/pum/ru/rut/Sh/ss/tko/Egfr/vn/w/Tm2/nAChRbeta2/up/ct/Rdl/Klp64D/mew/Mdr65/Gad1/Syn/pros/TkR99D/Su(H)/hop/14-3-3zeta/Orc3/Itpr/Galphaf/Lcch3/Nmdar1/cora/sinu/dia/Obp28a/Ms/Dop1R1/ine/Nos/cher/tim/drl/nAChRalpha3/Drip/nompA/vri/H15/Dyrk2/klg/mGluR/ppk/geko/CRMP/Clk/teq/CG3078/b6/Oamb/GluClalpha/Scgdelta/cry/Nipped-B/Eaat1/Inx7/Amph/CG4928/Rh50/nAChRalpha5/aru/NtR/CG5921/Pink1/CG5541/dpr18/Mec2/meng/Nha1/MICU1/CG8086/nAChRalpha6/dpr19/AstC/nrv3/Dscam1/Obp44a/KCNQ/CAP/ClC-b/Sans/dpr13/Obp56d/CG9313/Oatp58Dc/Eglp2/Eglp4/ITP/mri/dpr20/Chd64/DopEcR/CG10226/ppk26/TrpA1/Or71a/ClC-c/SecCl/Ir76b/SLC22A/CG14636/NPFR/Neurochondrin/mura/Desi/dpr5/Tk/ATP8B/dpr15/dpr9/subdued/Nrx-1/CG6972/Root/sosie/PQBP1/Obp99a/Obp99d/Cad99C/CG15537/Scamp/G9a/Alk/dpr1/robo3/Socs36E/Atx2/HDAC4/MCU/prt/dome/Proc/Obp83g/Obp83ef/Obp83cd/ClC-a/HEATR2/dpr17/Gba1b/wry/CG31935/dpr10/Rbfox1/Krn/be/dpr8/dpr11/dpr4/Nmdar2/dpr3/Dop2R/Unc-89/pain/Ctr1A/dpr12/Nox/NaCP60E/sif/mbm/kra/gish/DCX-EMAP/futsch/mnb/Pka-R1/brp/drd/qvr/Mp/2mit/rempA/Atg1/stj/pdm3/salm/Npc1b/Ank2/aPKC/mid/f/Shab/Tmhs/cac/spoon/dysc/orb2/CG43778/Ca-alpha1T/Mhc/Rgk1/pHCl-1/tn/rad/CrebB/rdgC/Fhos/rg/nAChRalpha4/SclA/SclB/Ggamma30A/Tmc/wupA/Phm/wg/Psn/para 241
GO:0006811 ion transport 231/3316 454/8880 0e+00 8.00e-07 8.0e-07 nAChRalpha1/nAChRbeta1/nAChRalpha2/Cam/eag/Ca-alpha1D/na/Atpalpha/sesB/Sh/Shaw/slo/tipE/w/nAChRbeta2/Syt1/Rdl/GluRIA/GluRIIA/Vha55/Itpr/Lcch3/Nmdar1/SerT/l(2)03659/ATPsynG/Elk/ine/Fer2LCH/nAChRalpha3/tty/nrv1/ATPsynCF6/ATPsynO/ppk/dve/Tsf1/CG3078/Picot/GluClalpha/homer/Mpcp2/Eaat2/Eaat1/CG4928/CG1688/Dic1/blot/GluRIID/JhI-21/VhaPPA1-1/Vha100-2/Vha100-1/Rh50/Nckx30C/CG4587/nAChRalpha5/CG15270/CG15279/NtR/SK/NAAT1/CG9657/CG10353/MFS10/Clic/CG11655/CG4239/CG9903/ppk23/CG7536/CG8034/CG8051/Usp2/MFS3/VGlut/OtopLb/CG10019/Nha1/MME1/MICU1/CG7781/CG14275/Mco1/nAChRalpha6/CG4995/Vha100-5/atilla/Oatp33Eb/spict/CG9336/CG9338/tadr/nrv3/CG8726/sand/Tmem63/CG8026/CG1698/KCNQ/Zip48C/CG8850/ClC-b/CG6329/VGAT/CG8468/DAT/CG9864/Mpcp1/Oatp58Dc/Tpc2/galene/CG5687/ZnT63C/Teh2/CG10483/Best2/ppk26/TrpA1/Ae2/CG6938/CG6928/Ncc69/Zip71B/CG7255/ClC-c/Papst2/SecCl/Oatp74D/CG5535/Sfxn2/SLC22A/CG1090/elm/CD98hc/Task7/Teh1/CG6723/sea/Task6/Zip88E/CG5404/Zip89B/ChT/CG7720/subdued/MFS9/KaiR1D/CG17119/Irk2/Mco3/gb/CG14507/Zip99C/CG9702/CG9717/ATPsynC/salt/Gat/CG11155/Best1/Vha36-3/olf186-F/MCU/Stim/GluRIIC/CG31030/ClC-a/GluRIIE/NKCC/Oatp26F/CG31668/CG31676/CG31935/Ndfip/CG32669/CG33181/ppk13/Nmdar2/dpr3/Mid1/pain/Ctr1A/EMRE/Rgk2/CG34396/Rgk3/NaCP60E/NKAIN/Orct2/CG42235/CG42237/CG42260/OtopLc/CG42340/Ca-beta/OtopLa/CG42496/qvr/CG42542/CG42594/stj/Syt12/Sytalpha/Sytbeta/CngA/Ca-Ma2d/Shab/CG43155/cac/Hk/Cngl/Ih/Vha68-2/Stacl/GluRIB/dysc/Slob/Ca-alpha1T/PsGEF/Rgk1/pHCl-1/Piezo/Trpm/Vha68-1/Neto/CngB/nAChRalpha4/Vha13/Psn/para 231
GO:0006812 cation transport 168/3316 319/8880 0e+00 6.00e-06 5.5e-06 nAChRalpha1/nAChRalpha2/Cam/eag/Ca-alpha1D/Atpalpha/Sh/Shaw/slo/tipE/w/Syt1/GluRIA/GluRIIA/Vha55/Itpr/SerT/ATPsynG/Elk/ine/Fer2LCH/nAChRalpha3/nrv1/ATPsynCF6/ATPsynO/ppk/dve/Tsf1/CG3078/homer/CG4928/CG1688/blot/GluRIID/JhI-21/VhaPPA1-1/Vha100-2/Vha100-1/Rh50/Nckx30C/CG4587/nAChRalpha5/CG15270/CG15279/SK/NAAT1/CG9657/CG4239/ppk23/Usp2/VGlut/OtopLb/Nha1/MME1/MICU1/CG7781/CG14275/Mco1/nAChRalpha6/CG4995/Vha100-5/atilla/spict/CG9336/CG9338/nrv3/sand/Tmem63/CG1698/KCNQ/Zip48C/CG8850/CG6329/VGAT/DAT/Tpc2/galene/CG5687/ZnT63C/Teh2/ppk26/TrpA1/Ncc69/Zip71B/CG7255/CG5535/Sfxn2/SLC22A/CG1090/elm/CD98hc/Task7/Teh1/CG6723/Task6/Zip88E/Zip89B/ChT/CG7720/subdued/KaiR1D/Irk2/Mco3/Zip99C/ATPsynC/salt/Gat/Vha36-3/olf186-F/MCU/Stim/GluRIIC/CG31030/GluRIIE/NKCC/CG31668/CG31676/Ndfip/CG32669/CG33181/ppk13/dpr3/Mid1/pain/Ctr1A/EMRE/Rgk2/CG34396/Rgk3/NaCP60E/NKAIN/Orct2/CG42235/CG42260/OtopLc/CG42340/Ca-beta/OtopLa/CG42496/qvr/CG42542/CG42594/stj/Syt12/Sytalpha/Sytbeta/CngA/Ca-Ma2d/Shab/CG43155/cac/Hk/Cngl/Ih/Vha68-2/Stacl/GluRIB/Ca-alpha1T/PsGEF/Rgk1/Piezo/Trpm/Vha68-1/Neto/CngB/Vha13/Psn/para 168
GO:0007610 behavior 238/3316 489/8880 1e-07 4.52e-05 4.2e-05 Abd-B/Appl/arm/b/ben/bsk/cin/clt/Ddc/disco/dnc/eag/for/Gld/Galphas/Gpdh1/Hn/kay/egh/dco/na/Atpalpha/Pkc98E/Poxn/pum/rut/Sh/slo/sm/ss/tipE/tko/Egfr/vn/w/oc/Syt1/Rdl/Gad1/5-HT7/Syn/pros/TkR99D/ogre/fru/Su(H)/hop/14-3-3zeta/ey/ple/Orc3/Itpr/Nmdar1/tutl/bol/Ms/Dop1R1/Mef2/EbpIII/Dh44/klu/Cdk5/Rac2/Rho1/cher/tim/Fer2LCH/dod/drl/Drip/vri/Pino/Dyrk2/klg/mGluR/ppk/lig/geko/Gpo1/Tsf1/Pka-R2/CRMP/Clk/aay/amon/Ac3/teq/Oamb/moody/cry/homer/tyf/Nipped-B/Eaat1/SPARC/Inx7/Cdk5alpha/Nep2/GABA-B-R2/qtc/Spn27A/aru/Surf1/SPR/Shmt/Pink1/BTBD9/CG5541/ppk23/CG6123/CCKLR-17D3/fipi/meng/Nlg2/MICU1/Dh31/AstC/sNPF/sand/trsn/CG13192/VGAT/Hr51/DAT/Obp56d/CG10433/ITP/MsR1/FMRFaR/DopEcR/ppk26/Pura/TrpA1/Tat/Prps/CrzR/Mip/Gem2/Gabat/sNPF-R/SLC22A/NPFR/mura/Tk/dpr9/CG8907/unc79/Ktl/SIFaR/Nrx-1/spas/Root/PQBP1/gb/unc80/jus/Obp99a/CDase/tara/krz/Scamp/G9a/Alk/dpr1/Socs36E/Atx2/HDAC4/rho-5/MCU/prt/dome/Ilp4/mtt/Tdc2/HEATR2/Gba1b/Trim9/Rbfox1/Krn/shep/be/Nmdar2/Dop2R/Mid1/pain/shakB/NaCP60E/Desat1/mbm/kra/gish/DCX-EMAP/futsch/mnb/Pka-R1/brp/Tdc1/qvr/Pdfr/2mit/verm/Ank2/aPKC/TfAP-2/Shab/dati/cac/Hk/Ih/Hmgcr/spoon/MsR2/dysc/orb2/Rgk1/Pde1c/Piezo/Nrg/tn/Neto/mbl/rad/CrebB/rg/lov/Dop1R2/wake/lncRNA:yar/unc-104/Tmc/Phm/Pal1/muc/wg/Sema5c/Psn/para/caz 238
Code
# Get the entrez IDs
up_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% up_genes, ]$entrezgene_id)

up_ekegg <- plotKEGG_dm(up_entrez, title = "KEGG, upregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "13 enrichments found"
Code
nrow(up_ekegg)
[1] 13
Code
write.csv(up_ekegg, file = paste0(output_dir, "O1O2KO_vs_ctrl_up_kegg.csv"))
Code
head(up_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme03430 dme03430 Mismatch repair - Drosophila melanogaster (fruit fly) 11/362 19/3488 0.0000005 0.0000511 0.0000456 Gnf1/PCNA/Pms2/spel1/tos/Pol31/CG8142/RfC3/Msh6/RfC4/Pol32 11
dme03420 dme03420 Nucleotide excision repair - Drosophila melanogaster (fruit fly) 18/362 51/3488 0.0000014 0.0000802 0.0000716 hay/mei-9/Xpac/Gnf1/PCNA/Pol31/Ercc1/CG8142/RfC3/Rpb11/Chrac-14/Rpb7/RfC4/mrn/Xpd/Rpb12/Cdk7/Pol32 18
dme03050 dme03050 Proteasome - Drosophila melanogaster (fruit fly) 18/362 53/3488 0.0000027 0.0001018 0.0000909 Prosbeta6/Prosalpha4/Prosbeta1/Prosalpha5/Prosbeta2/Prosalpha7/Prosbeta3/Rpn12/Prosbeta5/Prosbeta4/Pomp/Rpt3R/Rpt6R/CG30382/Prosalpha2/Prosbeta7/Prosalpha6/Prosalpha3 18
dme03440 dme03440 Homologous recombination - Drosophila melanogaster (fruit fly) 10/362 23/3488 0.0000425 0.0011894 0.0010619 mus101/spn-B/spn-D/mre11/Pol31/mms4/rad50/mus81/tefu/Pol32 10
dme03030 dme03030 DNA replication - Drosophila melanogaster (fruit fly) 12/362 34/3488 0.0000857 0.0019195 0.0017139 Gnf1/PCNA/DNApol-alpha73/Mcm5/Pol31/CG8142/CG13690/RfC3/Chrac-14/CG30105/RfC4/Pol32 12
dme03410 dme03410 Base excision repair - Drosophila melanogaster (fruit fly) 11/362 30/3488 0.0001142 0.0021046 0.0018791 Gnf1/PCNA/Pol31/CG8142/RfC3/CG9272/CG9601/Chrac-14/gkt/RfC4/Pol32 11
Code
# Get the entrez IDs
down_entrez <- na.omit(ensembl.genes[ensembl.genes$gene_id %in% down_genes, ]$entrezgene_id)

down_ekegg <- plotKEGG_dm(down_entrez, title = "KEGG, downregulated genes")
Running KEGG for organism = drosophila melanogaster

[1] "5 enrichments found"
Code
nrow(down_ekegg)
[1] 5
Code
write.csv(down_ekegg, file = paste0(output_dir, "O1O2KO_vs_ctrl_down_kegg.csv"))
Code
head(down_ekegg)
ID Description GeneRatio BgRatio pvalue p.adjust qvalue geneID Count
dme04080 dme04080 Neuroactive ligand-receptor interaction - Drosophila melanogaster (fruit fly) 36/963 67/3488 0.0000047 0.0005575 0.0005178 mAChR-A/alphaTry/Rdl/5-HT7/GluRIA/TkR99D/Lcch3/deltaTry/gammaTry/Nmdar1/etaTry/zetaTry/Dop1R1/Try29F/Lgr1/mGluR/AkhR/GABA-B-R2/CCKLR-17D3/GABA-B-R3/CCHa2-R/Dh44-R2/Dh44-R1/CrzR/mAChR-B/SIFaR/CCAP-R/lambdaTry/kappaTry/Dh31-R/Nmdar2/Dop2R/GABA-B-R1/5-HT1B/GluRIB/Dop1R2 36
dme04391 dme04391 Hippo signaling pathway - fly - Drosophila melanogaster (fruit fly) 30/963 56/3488 0.0000318 0.0018892 0.0017547 bsk/ft/hth/l(2)gl/dco/sd/tsh/vn/mts/ex/tws/14-3-3zeta/dia/Mad/14-3-3epsilon/par-6/wgn/lft/mats/fred/Krn/Patj/Pp2A-29B/hpo/aPKC/sdt/d/dally/wg/ds 30
dme04512 dme04512 ECM-receptor interaction - Drosophila melanogaster (fruit fly) 9/963 11/3488 0.0002833 0.0112383 0.0104381 Col4a1/LanA/vkg/Hml/CG3168/Tsp/Dg/wb/LanB2 9
dme04142 dme04142 Lysosome - Drosophila melanogaster (fruit fly) 49/963 123/3488 0.0018863 0.0561180 0.0521219 cm/Sap-r/DNaseII/g/Gal/or/CG13397/Rbcn-3A/Picot/AP-1mu/VhaPPA1-1/Vha100-2/Vha100-1/Tsp42Ed/Tsp42Ea/Ppt1/Gga/MFS10/Bace/LManIII/LManV/LManVI/Tsp29Fa/Lip4/Ppt2/Vha100-5/Lamp1/CG4847/CG13510/CG11459/Npc2b/CG17283/Sgsh/CG6656/CG17119/Npc2g/Npc2h/CG13516/Hexo2/CG18858/CG30269/Gba1b/CG33128/CG12163/Ect3/CG42565/CG42566/Npc1b/VhaAC45 49
dme00600 dme00600 Sphingolipid metabolism - Drosophila melanogaster (fruit fly) 16/963 31/3488 0.0038069 0.0906043 0.0841526 Sap-r/Gal/ifc/wun/CG33090/Sk1/nSMase/laza/CG11438/CG11437/CDase/Hexo2/bwa/Gba1b/Spt-I/Ect3 16

4 sessionInfo

Code
sessionInfo()
R version 4.2.2 (2022-10-31)
Platform: aarch64-apple-darwin20 (64-bit)
Running under: macOS Ventura 13.1

Matrix products: default
BLAS:   /Library/Frameworks/R.framework/Versions/4.2-arm64/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.2-arm64/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] grid      stats4    stats     graphics  grDevices utils     datasets 
[8] methods   base     

other attached packages:
 [1] org.Dm.eg.db_3.16.0         clusterProfiler_4.7.1.003  
 [3] gghighlight_0.4.0           eulerr_7.0.0               
 [5] forcats_1.0.0               stringr_1.5.0              
 [7] dplyr_1.1.3                 purrr_1.0.1                
 [9] readr_2.1.3                 tidyr_1.3.0                
[11] tibble_3.2.1                tidyverse_1.3.2            
[13] limma_3.54.1                ggrepel_0.9.3              
[15] geneplotter_1.76.0          annotate_1.76.0            
[17] XML_3.99-0.14               reshape_0.8.9              
[19] lattice_0.20-45             ggplot2_3.4.2              
[21] gridExtra_2.3               biomaRt_2.54.1             
[23] gplots_3.1.3                RColorBrewer_1.1-3         
[25] GenomicFeatures_1.50.4      AnnotationDbi_1.60.2       
[27] DESeq2_1.38.3               SummarizedExperiment_1.28.0
[29] Biobase_2.58.0              MatrixGenerics_1.10.0      
[31] matrixStats_1.0.0           GenomicRanges_1.50.2       
[33] GenomeInfoDb_1.34.9         IRanges_2.32.0             
[35] S4Vectors_0.36.2            BiocGenerics_0.44.0        

loaded via a namespace (and not attached):
  [1] utf8_1.2.3               tidyselect_1.2.0         RSQLite_2.3.1           
  [4] htmlwidgets_1.6.1        BiocParallel_1.32.6      scatterpie_0.2.1        
  [7] munsell_0.5.0            codetools_0.2-19         withr_2.5.0             
 [10] colorspace_2.1-0         GOSemSim_2.24.0          filelock_1.0.2          
 [13] knitr_1.42               rstudioapi_0.14          DOSE_3.24.2             
 [16] labeling_0.4.2           bbmle_1.0.25             GenomeInfoDbData_1.2.9  
 [19] polyclip_1.10-4          bit64_4.0.5              farver_2.1.1            
 [22] downloader_0.4           coda_0.19-4              vctrs_0.6.3             
 [25] treeio_1.22.0            generics_0.1.3           gson_0.1.0              
 [28] xfun_0.37                timechange_0.2.0         BiocFileCache_2.6.1     
 [31] R6_2.5.1                 apeglm_1.20.0            graphlayouts_1.0.0      
 [34] locfit_1.5-9.7           bitops_1.0-7             cachem_1.0.8            
 [37] fgsea_1.24.0             gridGraphics_0.5-1       DelayedArray_0.24.0     
 [40] BiocIO_1.8.0             scales_1.2.1             ggraph_2.1.0            
 [43] enrichplot_1.18.4        googlesheets4_1.0.1      gtable_0.3.3            
 [46] tidygraph_1.2.3          rlang_1.1.1              splines_4.2.2           
 [49] rtracklayer_1.58.0       lazyeval_0.2.2           gargle_1.3.0            
 [52] broom_1.0.3              yaml_2.3.7               reshape2_1.4.4          
 [55] modelr_0.1.10            backports_1.4.1          qvalue_2.30.0           
 [58] tools_4.2.2              ggplotify_0.1.1          Rcpp_1.0.11             
 [61] plyr_1.8.8               progress_1.2.2           zlibbioc_1.44.0         
 [64] RCurl_1.98-1.12          prettyunits_1.1.1        viridis_0.6.3           
 [67] cowplot_1.1.1            haven_2.5.1              fs_1.6.2                
 [70] magrittr_2.0.3           data.table_1.14.8        openxlsx_4.2.5.2        
 [73] reprex_2.0.2             mvtnorm_1.1-3            googledrive_2.0.0       
 [76] hms_1.1.3                patchwork_1.1.2          evaluate_0.20           
 [79] xtable_1.8-4             HDO.db_0.99.1            emdbook_1.3.12          
 [82] readxl_1.4.1             bdsmatrix_1.3-6          compiler_4.2.2          
 [85] KernSmooth_2.23-20       crayon_1.5.2             shadowtext_0.1.2        
 [88] htmltools_0.5.4          ggfun_0.1.1              tzdb_0.3.0              
 [91] aplot_0.1.10             lubridate_1.9.1          DBI_1.1.3               
 [94] tweenr_2.0.2             dbplyr_2.3.3             MASS_7.3-58.2           
 [97] rappdirs_0.3.3           Matrix_1.5-4.1           cli_3.6.1               
[100] parallel_4.2.2           igraph_1.5.0             pkgconfig_2.0.3         
[103] GenomicAlignments_1.34.1 numDeriv_2016.8-1.1      xml2_1.3.5              
[106] ggtree_3.6.2             XVector_0.38.0           rvest_1.0.3             
[109] yulab.utils_0.0.6        digest_0.6.33            Biostrings_2.66.0       
[112] rmarkdown_2.20           cellranger_1.1.0         fastmatch_1.1-3         
[115] tidytree_0.4.4           restfulr_0.0.15          curl_5.0.1              
[118] Rsamtools_2.14.0         gtools_3.9.4             rjson_0.2.21            
[121] lifecycle_1.0.3          nlme_3.1-162             jsonlite_1.8.7          
[124] viridisLite_0.4.2        fansi_1.0.4              pillar_1.9.0            
[127] KEGGREST_1.38.0          fastmap_1.1.1            httr_1.4.6              
[130] GO.db_3.16.0             glue_1.6.2               zip_2.2.2               
[133] png_0.1-8                bit_4.0.5                ggforce_0.4.1           
[136] stringi_1.7.12           blob_1.2.4               caTools_1.18.2          
[139] memoise_2.0.1            ape_5.7-1